{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": 1,
   "id": "a684a3c3",
   "metadata": {},
   "outputs": [
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "Matplotlib is building the font cache; this may take a moment.\n"
     ]
    },
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "1.8517857142857121\n",
      "2.389642857142869\n",
      "2.4678571428571336\n"
     ]
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAW4AAAEGCAYAAABFBX+4AAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjUuMSwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/YYfK9AAAACXBIWXMAAAsTAAALEwEAmpwYAAAyjklEQVR4nO3dd3hU1dbH8e+CBAhFqaIBBRVFmoCGYggQqVIUghSVYqHIBeSCoIhw9b6KyDWIQSNNUBQQRUQEAyI1EnroSJEqXSLSCYQk6/1jJpGeAEnOTLI+z8NzJzOnrJnH+5s9++x9tqgqxhhjvEc2pwswxhhzYyy4jTHGy1hwG2OMl7HgNsYYL2PBbYwxXsaC2xhjvIwFt8myRCRPOh77tvQ6tjEW3MbjiMgAEXn2Brb3ERGfqzyfQ0TkOrv+LiKPi0hXEfG/yv7zRaTqRX+3EJGG7seviUil6xx7hoj0vk7NPUQkl4gsEJEK7uPdJiKjRKTWNfZpIyL3XuecJou44j92Y9KKiEwCygCJl72UDfhVVXtdZZ9OwDvAShHpftFLX6vqCPc2QUAj4ALgC/wBdBCReKAqsAU4BeQEnnG/fjWn3ceoALQAGlz2+gXgjDv8Q4HvgPEiUsF9nlXXeN+9gCNACxGJVNU1V9nMBxgIxAN5gGdUNVRE6gIfXKPeJsA+YLeIHAc2XfRaNiBeVa8a+iZzseA26UZV26Z2WxEpCLwP5AJqqOpy9/MBwHhg0kWb7wSW4wr4zsA6VR3j3v4voJ6qHrvo2NkBUdX4y057BlfgvQ70FJHcqnpWRJoC7wIlgG+AN3EFa18R2Qt0ABLc/y5/Hy2BZ4E6QGFgpoh0U9Woi7bJBRwEFCgABAO/ulvaeVV1l3u7nKp6XkTyurfLCRQUkSLAblUNuuiYeYGfUviYTSZhwW08RRtgG/AREOVura8FPgdaq+qJpA1V9ZCIbATOqGo0gIhEXXSsmSKSEwhV1SlAY+AjEUnAFZZJigFf42p5ZwPaA2VV9SfgJxH5CeilqjtEJNS9TzfgKFCfi7oa3V8OrwItgSaqegZXa/1pXN0mPwP/U9XDQA6gFq5WewCuXwh/A/8CsolINHAPcMrdHVMZ6OHeNj9QGijmfs+Cq8V+Fjh5Ix+48V4W3CZV3H3I2YE4vcYNbtxdCjlw/fyfj6ub5AJwO3AC139vOXCFTHZgiqr+G0BVR150nPbAOmAP0E9VN6aixPtU1V9EfFQ13t1dkcd97JnAzKvUOwLYlNQFcw2PiEjupD9Udad7X1/cLW4RuR9XN8p2XF8+X4tIUveQAL/g+sLYLiKPubfZheuL4xxwHDgGBALhqvquiHwDfKiqp4BfcbXIxwPjcAX0s6oaJCIFgAhVDUzFZ2QyCQtuk1rtgC8Arn+9D4AQVQ12b5sdV8u4uLsL4hlVbXf5Du7wewRXP3Mg8B6u4G8HzEhFfUlBOUFEprgfJ3/BiEi4+7g5cQXuWVwhGHCVWooCbd2vHQOWXfRaRVVdzz9fRuD6gumvqnPc3RgXn1uARFU9KiLvq2qM+wLjX7i6aN4B3gaa4+qvf9S93z24uoSSzuuDq+tlHK6+8QfcLfNsFz0W4BNVHZ+Kz8t4MQtuk1rTgCggjisvNibJhqtFffii5wrialGmJARXWE3B1ac8GigHxLi7BAq5j19WVZP7lkWkBK6uiyS5cAVpiYsPrqo93NuvA5qr6h4RKQm8KSJy2a+IB4DcwFbg/1R1n4ggIk2A90TkUVzdFbvcx04A5ojIQFwBfAUR+VhVv7roM+kN/J/7Pc3E9cskEFgkIvmBfKr6t3vfQFxfmrlw9a+vwdWvXzupb1tVg92/eLJf8xM2mYYFt0kVVT3JzfWh1gVWpuL4Q91dEh2Az4BvgfdUdbU7kBYAH6tqgntUx4tAFfd2T/JPf/M9uEZb1MbVAr3eOfeIyAlcFwcXXvR8FK5+9uqQ3Nq9CxgENAQeB4oA/ri6R5L4A4NUdfrF5xGR/+L64kk6/mr3yJgooLr7mI1VNVZEpgE/ut9vkq3AU0B/4DxwN7BfRCKA24AKIrII2KWqL13vPZvMwcZxm3Tj7gf+NzDxGq8/cNk46/ZAeaAmrouSE919wuOAtar6g3u7wri6MO5T1eqqGoOrxdsQ2A30xdUVsucapZUXkZ7uluxg4P/EpZOIlLvK9nfguoDY1P34U2AAMNV9ETTJtX6JwKUXRcEV8quBn3H9uljtfn4arguXyX3yqvq3qm67aN/q7u2LqmpNVc3v7pp6+DrnN5mIBbdJF+KalfgFcMQ9sgNcwXaXuCbM+AIRuEZnAKCqo1W1h6ruwzX2+htcLeGdqvrqRdstVNV3VfWA+1y+uIYHDgF64up7fhtX6xQRySkiQ9wjO/Li6s44iau1/B2ulvlHuFrxtd37vIBrNMcZVT2Ia3hfZ1xB+y9VHYZr1Mu/Ln7bwFARib74H9DlKh/RDly/FnLi+hIoJyIt3M/1BcaISLXL9smG68LvE9jQvyzNukpMmnN3eSwFfsfVik6yBvDD1ZWhuLoK5rn3qYertV0CqAgUxRWqD6rq3uucqzquroUfgNruLp0RIrIZ+FBEarjHQp/CdcEz6vJRMSISgmtYYCn3ccA1PK+Nqh5y/30PUBKopqr73c/1xTUWPIkv0PcaXSW+l30+EbhGzrRQ1cMi0hbXUMNnVPV3EVkNvCsizd1jy7vh6r8+A6h7m6QhgUkubv2bTExs6TKTHkSkiLsLI7XbB+IaUbITVzfAlmsNO7zKvne7W+mXP59NVa/XfeHVRMRPVWOdrsNkPAtuY4zxMtbHbYwxXibd+7gLFy6sJUuWTO/TGGNMprJ69eq/VLXI1V5L9+AuWbIk0dHR6X0aY4zJVETkWne1tK4SY4zxNhbcxhjjZSy4jTHGyzgyAefChQvs37+fc+fOOXF6r5IrVy6KFy+Or69vyhsbY7IER4J7//795MuXj5IlS6bmFqFZlqpy9OhR9u/fz7332lKDxhgXR4L73LlzFtqpICIUKlSImJhUT0A0xniA6WsPEDpnGwePx+Kf34/XGpameeViaXZ8x+5VYqGdOvY5GeNdpq89QP9pG4m94Lpt/IHjsfSf5lrEKa3C2y5OGmNMGgqdsy05tJPEXkggdM62a+xx4yy4b8CJEydI7b1dEhISOH/+fDpXZIzxNAePu+77dXLVdGJ3Rl/xfFpIdXCLyAgRedL9eJyILHUv1eS19u3bR0hISPLfn376KQcPHmTNmjVMnHjlvf+feeYZhg4det1jtmzZks2bN7Nw4UJ69ep1xesxMTG0aNEi+e8OHTqwYsUKvvnmm5t/I8YYj+Gf349TayI4tmAsZ7ZEXvJ8WklVcItITeBOVZ3pvtl7dveq0v4i8kCaVZOB4uLi8PHxIXfu3Fy4cIE+ffpQvHhxBg8ezLFjx9ixY8cl24eGhlKrVi22bt3Kzz//fM3jJh3Tx8eH7NmvXP5vxIgR9O7dm927dwPg6+tLgQIF+OOPPzhw4EDavkljTIYLuLCRv+eOxO/+KhRq1BMAP9/svNawdJqdI8WLk+7VRT4DZolIM1zr8yWtaLIACOLSdfduSK9evVi3bt3N7n5VlSpVIiws7LrbTJ8+nREjRrB9+3aGDh3Kpk2bCA0NZciQITRt2pRs2f75Tvvwww/ZuXMno0aNIi4ujrZt27JhwwZ69+6Nr68vBw8e5KGHHuKRRx5h8+bNPPPMM5w7d46YmBg2bdrEunXrWLNmDefPn+fo0aPkzJmTF154gf/85z9s376d/v37kz9/ft59911GjRqVpp+FMSbjfP/994T/91UqVKlBnicHcvhMgmOjSjoAm4EPgFeA7rjWAATX8k+lLt9BRLrgXq7pnnvuSZNC01rr1q2pU6cOffv2pV+/fixZsoRs2bIRFRXFwoULyZYtGzt37qR79+7ky5eP8uXLM2jQIADKlSvHzp07qVSpEqNGjeKhhx6ibNmyLFq0iGeeeYYhQ4awZ88epk6dSnh4OMHBweTIkYPPP/+cpUuX8t1337Fo0SLGjh3LnXfeSc+ePalRowbr169n7969HvuZGWOubdasWTz77LNUr16dOXN+Jm/evOl2rtQEd2VgjHt5pYlAIK7lp8C1ft8V3S2qOgYYAxAQEHDdq3kptYzT04YNG4iMjGTEiBHJzyUmJnL27Fly5cpF3rx5GTBgAOXLl2fPnj3JQ/MSExN54IEHOHPmDPny5SM2NpbNmzcTHBx8zRY3wKBBg2jRogVz586lVKlSlCxZkpUrVzJgwAD8/f25++67eeklW6TbGG+zaNEinn76acqXL09ERES6hjakLrh3APe5HwfgWncvCNfirBWBtBvjkoEmTJjAsmXLqFWrFj169ODnn39m9uzZrFixgttvv517772XokWLUqRIEcqXL4+/v/8l+8fGxrJkyZLkxym1uAG2bdtG69atad++Pd9++y2lS5dm4MCBREREULVqVRYvXsxDDz2U0R+FMeYWLF++nKZNm3Lffffxyy+/kD9//nQ/Z2qCexzwuYg8g2vB02Bghoj4A42A6ulXXvpp3749jRo1om/fvoBrhMno0aP5+uuvCQ4OplOnTgBky5aNO++8k3nz5l2yf/XqN/62CxQowNChQ7nvvvvIly8f9957LwkJCXTv3p1t27bx/fff3/obM8ZkmHXr1tGoUSPuvPNO5s6dS+HChTPkvCkGt6qeAlpd/JyIBAP1gQ9U9US6VJYBksZkJyQk8MQTT/D666/TvXt3WrVqxaRJk6hduzalS5dmw4YN1KtX75J9jx07dslxUtNVsmbNGiZNmkRCQgIFCxZkyJAhDBgwgJIlS/Lggw+yfft2tm/fTsuWLTPsMzDG3JytW7fSoEED8ubNy/z586/4VZ6uVDVd/z366KN6uc2bN1/xXEZLSEjQxx9/XEeOHKmqqtu2bdOyZctqeHi4qqquXLlSP//8c01ISNDg4OAr9n/kkUeSH+/fv1+rVat2zXPVrl1bd+3apUeOHNE9e/aoqurevXu1QoUKOnPmTFVV/e9//6uFChXSRYsWXbG/J3xexph/7Ny5U/39/bVo0aK6bdu2dDkHEK3XyNV0X+U9ICBAL1+6bMuWLZQpUyZdz5saqpp8wVFVOXz4MHfddZcj5wfXzbdy5cp1xXae8nkZY1x3N61VqxbHjx8nMjKSChUqpMt5RGS1qgZc7TXHbjLlCS4OTRHJ0NC+/PzAVUPbGOM5jhw5Qr169fjrr7+YP39+uoV2SrJ0cBtjTGodO3aMBg0asHfvXn7++WeqVKniWC0W3MYYk4JTp07RqFEjtmzZwowZM6hVq5aj9VhwG2PMdcTGxvLUU08RHR3N1KlTadiwodMlZd3gTkxMBEi+J8n//vc/SpUqxdNPPw241sXMli0b3333HSdOnGD//v0EBwdTt25dOnbsyCuvvEKlSpWSj/fWW29Rp04d5s6dS758+ejevTstW7Zk1qxZl9xsKj4+ntatWzNt2jSee+45Dh48mPxa3rx5+emnnzLg3RtjUiMuLo6nn36ayMhIJkyYQPPmzZ0uCfCS4E6PZYA2btxInz598PHx4ciRIxQpUoQFCxbQp08fSpUqhaoyaNAgwsPDSUhI4MyZMyxZsoQKFSqwdOlSxo4dmzw05+zZs9x2220sXbqUI0eO8Oeff7Jnzx7y5MlD9uzZk78kRAQfHx/i4uIA15fDokWLkmuy8dvGeI74+Hiee+45Zs+ezZgxY2jbtq3TJSXz+OBOr2WAKlasSMeOHRERcubMyVtvvcX58+dp0qQJffv2JTw8nIcffpju3bsnt7gfe+wxoqOjiYuL4/HHH2ffvn2EhoZSrVo1jh49Snh4OJUqVSIwMJBPP/2UHTt2UKtWLXbs2MH06dPZt28fY8aMYd26dbRs2ZJ9+/YRHBxMQkICIsL+/fvT5DMzxtyaxMREOnbsyPfff8+wYcPo3Lmz0yVdwuOD+3rLAN1qqzshIYE9e/awceNG3nzzTe6//34+++wzZs+ezZkzZ1ixYgWHDx9GRChcuDBbtmzhiy++oGnTpjz33HNERkbSokULDh8+zJ49e+jTpw9bt27l8OHDbNiwgffee49SpUoxevRoqlatStWqVXn66adp2rQp3bp1Y8yYMXzzzTdMnjwZgGefffaW3o8x5tapKj169OCrr77inXfeoXfv3k6XdAWPD+5rLfeTFssAVa9enfj4eIoXL86XX35J7dq1qVGjBgDt2rUjKCiIwYMHEx8fT2JiIj4+PgwfPpzx48ezY8cO7r77bsD1k+rtt99mwoQJ9OnTB19fX958801Wr15Njhw5uO+++5LPuWjRIhYvXkxISAhbt26lXr16HDlyBIDRo0fz3HPP0aVLl1t+b8aYG6eq9OvXj5EjR/L6668zcKBnLvLl8cHtn9+PA1cJ6bRYBqhLly4kJiaiqqxdu5bY2NjkSTHx8fH88MMP5MqVi3r16hEfH8+8efOoV68eP/30E+PHj+err74C4I8//mDgwIFs376d9evXs27dOnbv3p28RFmjRo0A+Oqrr/j++++pWbMmHTt2ZOHChUycOJGpU6cC1sdtjNMGDRpEaGgo3bp1Y8iQIVdMkvMUHr9Y8GsNS+Pne+kSYGm1DFDu3LmZP38+DzzwABEREURGRhIYGMisWbMoX7588mK/06dPZ8aMGcn7Va1alV27dlGgQAEAatSoQevWralbty7t2rWjbdu2+Pr68sgjj/Djjz8mD9Rv1aoV06dPB2Dv3r0UKVKE+++/n/DwcMLDwylWLO1WyDDG3Jhhw4bx1ltv0aFDBz755BOPDW3wghZ3Uj92Wo8qSRIZGcmcOXPIlSsXM2bM4PTp08l9ziKCqrJo0SJOnz5Ny5YtiYmJITw8nK5du9KpUyfGjx+Pj48PXbt2JU+ePAwaNIgSJUqwe/dufvvtN3LkyMGaNWsICAjAz88v+bhRUVHUqVOHQ4cOJS8U/Mwzz6TJezLG3JgxY8bQp08fWrZsybhx4y5ZutAjXevuU2n1z1PvDqiqWrduXU1ISNDExMRLnh8yZIhWqFBBz549q3Xr1lVV1Z49e2qfPn20cePGum7dOlVVffPNN3XKlCn6999/a5s2bbR///566tQpXbFihdauXVtXrVqle/fu1apVq+qOHTtUVfX111/Xjz/+WLt166ZxcXF67733au3atbV27dp65513XrVOT/m8jMmMJk6cqCKijRs31vPnzztdTjLs7oA3Ri+7a9/FEhMTr/ptHB8fj4+PT/L+iYmJyRNvrne8Cxcu4Ovre916PP3zMsZb/fDDD7Rq1YqaNWsya9as5F/FnuB6dwd07PdAen9h3Irr9W1d6ydUUmgn7X/xbMnrHS+l0Pbkz8kYbzZnzhzatGlDlSpVmDFjhkeFdkocCe5cuXJx9OhRC6UUqCpHjx61270ak8Z+/fVXQkJCKFeuHLNmzSJfvnxOl3RDHLk4Wbx4cfbv309MTIwTp/cquXLlonjx4k6XYUymsXLlSpo2bUqJEiWYM2dO8ugwb+JIcPv6+nLvvfc6cWpjTBa2YcMGnnjiCQoXLsy8efO44447nC7ppnj4mBdjjEkbv//+O/Xr10+ev+HN8yYsuI0xmd6ePXuoW7cuqsq8efO8/he/x0/AMcaYW3Hw4EHq1q3L6dOnWbRoEQ899JDTJd0yC25jTKYVExND/fr1OXLkCPPmzaNixYpOl5QmLLiNMZnS8ePHadiwIbt27WL27NlUq1bN6ZLSjAW3MSbTOX36NI0bN2bTpk38+OOPBAcHO11SmrpucIuID7DL/Q/gFaAl0BhYoao90rc8Y4y5MefOnaNZs2asWLGCKVOmJN9WOTNJaVTJw8BkVQ1W1WAgJxAEVAX2i0i9dK7PGGNSLS4ujpYtW7JgwQLGjx+fvPh3ZpNScFcHQkQkSkQmAXWA7913rpoH1LzaTiLSRUSiRSTaZkcaYzJCQkIC7du3JyIigpEjR9K+fXunS0o3KQX3KqC2qgYBxwE/4ID7tZNA0avtpKpjVDVAVQOKFCmSVrUaY8xVJSYm0qlTJ6ZMmUJoaChdu3Z1uqR0ldLFyQ2qet79eCuQA1d4A+TFJvAYYxymqvz73/9m/PjxvP322/Tt29fpktJdSsE7QUQqikh2IATIg6uPG6AisCcdazPGmBQNGDCA8PBw+vTpw9tvv+10ORkipRb3O8DXgAAzgEHAYhEZDjzh/meMMY4YPHgw77//Pi+//DKhoaEevU5kWrpucKvqJlwjS5K5R5I0AYar6u50rM0YY65p+PDhDBgwgHbt2jFixIgsE9pwExNwVDUWmJoOtRhjTKqMGzeOXr16ERISwhdffOH5i/umsaz1bo0xXm/y5Ml07tyZJ554gsmTJ1+ybGBWYcFtjPEaM2bMoH379tSsWZPvv/+enDlzOl2SIyy4jTFeYd68ebRq1YpHH32Un376idy5cztdkmMsuI0xHi8qKopmzZrx0EMPMXv2bK9b3DetWXAbYzza6tWradKkCcWLF+eXX36hYMGCTpfkOAtuY4zH2rRpEw0aNKBgwYLMnz+fokWvepeNLMeC2xjjkbZv3079+vXJmTMn8+bNo3jx4k6X5DGy3jgaY4zH27t3L3Xr1iU+Pp7IyEjuv/9+p0vyKBbcxhiPcujQIerWrcvJkydZuHAhZcuWdbokj2PBbYzxGEePHqV+/focOnSIuXPnUrlyZadL8kgW3MYYj3DixAkaNmzIjh07mDVrFo899pjTJXksC25jjOPOnDlD06ZNWb9+PdOnT6dOnTpOl+TRLLiNMY46d+4cISEhLF26lG+++YYmTZo4XZLHs+A2xjjmwoULtGnThrlz5zJ+/HhatWrldElewcZxG2MckZCQwPPPP8+MGTMIDw/n+eefd7okr2HBbYzJcImJibz88stMnjyZ//3vf3Tv3t3pkryKBbcxJkOpKq+++irjxo1j4MCBvP76606X5HUsuI0xGeqtt95i+PDh9OrVi3feecfpcrySBbcxJsMMGTKEQYMG0alTJ4YNG5al1olMSxbcxpgMER4eTv/+/XnuuecYNWqUhfYtsOA2xqS7sWPH8sorr9CsWTPGjx9P9uzZnS7Jq1lwG2PSTWJiIi079qRz58743fsofwZ0JWLTEafL8no2AccYky7Onj1LvWatWTYvgrwVn6Bg/a4cOp1A/2kbAWheuZjDFXova3EbY9LcoUOHCA4OZtn8WRR4vCMFG3ZHsrvaibEXEgids83hCr1bqoJbRIqKyFr343EislREBqZvacYYb7R+/XqqVavG5s2buSNkILdVDbniQuTB47EOVZc5pLbFPRTwE5EWQHZVDQT8ReSB9CvNGONtZs6cSY0aNVBVoqKiKFUl+Krb+ef3y9jCMpkUg1tE6gBngMNAMDDF/dICICjdKjPGeA1V5aOPPqJZs2aUKVOGlStXUqlSJV5rWBo/30tHkPj5Zue1hqUdqjRzuG5wi0gO4C3gDfdTeYAD7scngasuuSwiXUQkWkSiY2Ji0qpWY4wHunDhAl27duXVV18lJCSEyMhI7rrrLsB1AfL9FhUolt8PAYrl9+P9FhXswuQtSmlUyRvAp6p63N1HdRpI+o2Tl2sEv6qOAcYABAQEaNqUaozxNMePH6dVq1bMmzeP/v37M2jQILJluzQWmlcuZkGdxlIK7npAHRHpDlQC7gH2AcuBioBdGjYmi9q5cydNmzZl586dfPHFF7zwwgtOl5RlXDe4VbVW0mMRWQQ8BSwWEX+gEVA9XaszxnikxYsXExISgqoyb948atWqlfJOJs2kehy3qgar6klcFyiXA4+r6on0KswY45m++uor6tatS6FChVixYoWFtgNueAKOqh5T1Smqejg9CjLGeKbExEQGDhzI888/T82aNVm+fDmlSpVyuqwsyaa8G2NSdPbsWZ5//nmmTp1K586d+fTTT/H19XW6rCzLgtsYc12HDh2iWbNmREdHM3ToUF599VW7JavDLLiNMde0fv16nnzySf7++2+mT5/OU0895XRJBrvJlDHmGn766SeCgoJITExk8eLFFtoexILbGHMJVSUsLIxmzZpRunRpVq5cSeXKlZ0uy1zEgtsYk+zChQv861//onfv3jRv3pxff/0Vf39/p8syl7HgNsYArunrTZo0YfTo0bzxxht899135M6d2+myzFXYxUljzCXT1z///HNefPFFp0sy12HBbUwWFxUVRfPmzVFV5s6dS+3atZ0uyaTAukqMycImTJiQPH19+fLlFtpewoLbmCwoafp6hw4dqFGjBsuWLeOBB2xBK29hXSXGZDGxsbE8//zzfPfdd3Tq1IkRI0bY9HUvY8FtTBZy+PBhmjVrxqpVq2z6uhez4DYmi9iwYQNNmzbl6NGj/PDDDzRr1szpksxNsj5uY7KAiIgIatSoQWJiIlFRURbaXs6C25hMTFUZPnw4Tz31lE1fz0QsuI3JpC5cuEC3bt3o1asXzZs3JzIy0qavZxIW3MZkQknT10eNGkW/fv347rvvyJMnj9NlmTRiFyeNyWR27dpF06ZN2bFjh01fz6QsuI3JRKKioggJCSExMdGmr2di1lViTCYxceJE6tatS8GCBW36eiZnwW2Ml0tMTOQ///kP7du3t+nrWYR1lRjjxWJjY3nhhReYMmUKHTt2ZMSIEeTIkcPpskw6s+A2xktdPH09NDSUPn362PT1LCJVwS0iBYFHgbWq+lf6lmSMScmGDRt48skn+euvv5g2bRrNmzd3uiSTgVLs4xaRu4AIoCqwUESKiMg4EVkqIgPTvUJjzCWSpq/Hx8ezePFiC+0sKDUXJ8sBvVX1PWAOUAfIrqqBgL+I2FUQYzLAxdPXH3zwQVauXMkjjzzidFnGASl2lajqPAARqYWr1V0QmOJ+eQEQBGy/eB8R6QJ0AbjnnnvSsFxjsqb4+Hh69uzJyJEjCQkJYcKECTYTMgtL1XBAcV3xaANcAAQ44H7pJFD08u1VdYyqBqhqQJEiRdKqVmOypKTp6yNHjqRfv35MnTrVQjuLS9XFSVVVoLuIvAu0BPzcL+XFxoIbk26Spq9v376dcePG8dJLLzldkvEAKQa3iPQDDqnqV0B+YAiu7pHlQEVgW3oWaExWtWTJEpo3b05CQgJz584lODjY6ZKMh0hNa3kM0F5EfgWyA9Pdfw8DWuMacWKMSUOTJk2iTp06FChQgBUrVlhom0ukGNyqekxV66tqLVXtpqongGBcLe7H3X8bY9JAYmIib731Fu3atSMwMJDly5fb9HVzhZuaOamqx/hnZIkxJg1cPH39pZdeYuTIkTZ93VyVTXk3JgNNX3uA0DnbOHg8Fv/8frzWsDTNKxfj8OHDNG/enJUrV/LBBx/Qt29fm75ursmC25gMMn3tAfpP20jshQQADhyPpf+0jfyxfQvDXuto09dNqtlQPmMySOicbcmhneTvrct5tf1TNn3d3BBrcRuTQQ4ej01+rKqcWj2TYwvGkuOOe1m5MpJixYo5WJ3xJhbcxmQQ//x+HDgeiyYm8Pe8MZxeG4HfA9V5uN1AC21zQ6yrxJgM8lrD0mQ7vp8/J/Xj9NoIbqv2NPe0/g9vPFnJ6dKMl7EWtzEZIC4ujvUzxrJv3Hvgm4vCT77Gg4FPJI8qMeZGWHAbk86WLVtG586d+e2333juuef46KOPuOOOO5wuy3gx6yoxJp2cOnWKnj17UqNGDU6ePElERASTJk2y0Da3zILbmHQwe/ZsypcvT3h4OD169OC3336jcePGTpdlMgkLbmPSUExMDO3ataNx48bkyZOHqKgoPv74Y/Lly+d0aSYTseA2Jg2oKhMnTqRMmTJMmTKFt99+m7Vr1xIYGOh0aSYTsouTxtyiP/74g65du/Lzzz9TvXp1xo4dS7ly5Zwuy2Ri1uI25iYlJCTw8ccfU65cORYvXszHH39MVFSUhbZJd9biNuYm/Pbbb3Ts2JEVK1bQqFEjRo0aZQtjmwxjLW5jbsD58+d5++23qVy5Mjt37mTSpElERERYaJsMZS1uY1Jp6dKldOrUiS1bttCuXTs++ugjChcu7HRZJguyFrcxKTh16hQ9evQgKCiIM2fOMHv2bCZMmGChbRxjwW3MdURERFC2bFlGjBhBz549+e2333jiiSecLstkcRbcxlzFkSNHePbZZ2natCm33347S5cuJSwsjLx58zpdmjEW3MZcTFX56quvKFOmDNOmTeOdd95hzZo1VK9e3enSjElmFyeNcdu9ezddu3bll19+ITAwkM8++4yyZcs6XZYxV7AWt8nyEhIS+OijjyhfvjxLly4lPDycxYsXW2gbj2UtbpOlbdiwgU6dOrFq1SqaNGnCyJEjufvuu50uy5jrsha3yZLOnTvHwIEDefTRR9mzZw+TJ09m5syZFtrGK6TY4haR24Fv3NueBtoAI4EywCxVHZSuFRqTxhYvXkznzp3Ztm0bHTp0YNiwYRQqVMjpsoxJtdS0uNsCw1S1PnAYeAbIrqqBgL+IPJCeBRqTVk6ePEm3bt2oVasW58+fZ86cOXz55ZcW2sbrpBjcqjpCVee6/ywCtAOmuP9eAARdvo+IdBGRaBGJjomJSbNijblZM2fOpGzZsowePZrevXuzadMmGjRo4HRZxtyUVPdxi8hjQAFgH3DA/fRJoOjl26rqGFUNUNWAIkWKpEmhxtyMP//8kzZt2vDUU09RsGBBli1bxrBhw8iTJ4/TpRlz01IV3CJSEPgEeAlXP7ef+6W8qT2GMRlJVRk/fjxlypRh+vTpDBo0iOjoaKpWrep0acbcshRDV0Ry4Ooa6a+qfwCr+ad7pCKwJ92qM+Ym7Nq1iwYNGvDiiy9Srlw51q9fz4ABA8iRI4fTpRmTJlLTWu4IPAoMEJFFgADtRWQY0BqISL/yjEm9+Ph4PvzwQ8qXL8+KFSsYOXIkkZGRPPTQQ06XZkyaSnE4oKqOxDX8L5mIzADqAx+o6ol0qs2YVFu/fj2dOnUiOjqaJ598khEjRlC8eHGnyzImXdxU/7SqHlPVKap6OK0LMuZGnDt3jgEDBhAQEMDevXv59ttv+fHHHy20TaZmU96N14qMjKRLly78/vvvvPDCC3z44YcULFjQ6bKMSXc2IsR4nRMnTvDyyy8THBzMhQsXmDt3Ll988YWFtskyLLiNV5k+fTply5Zl7Nix9O3bl40bN1KvXj2nyzImQ1lwG69w+PBhWrVqRUhICIULF2bFihWEhobaRBqTJVlwG4+mqowbN44yZcowc+ZMBg8eTHR0NAEBAU6XZoxj7OKk8Vg7duzg5ZdfZsGCBdSqVYvPPvuMBx980OmyjHGctbiNx4mPjyc0NJQKFSoQHR3N6NGjWbhwoYW2MW7W4jYeZe3atXTq1Ik1a9bQrFkzPv30U4oVK+Z0WcZ4FGtxG48QGxvLG2+8QZUqVThw4ABTp07lhx9+sNA25iqsxW0ct2jRIjp37syOHTvo2LEjoaGhFChQwOmyjPFY1uI2jlBVFi5cSLNmzXj88cdJTExk/vz5jB071kLbmBRYi9tkqHPnzjF58mTCwsLYsGEDPnluJ3+N5yjSoC0nC5R2ujxjvIIFt8kQhw8fZuTIkYwcOZKYmBhKlHqIO5v2IkfpWohPDg6fVfpP2whA88rWr23M9VhXiUlX69at44UXXqBEiRK8++67VK9enfnz5+P/Ujg5y9VDfP5Z3CD2QgKhc7Y5WK0x3sFa3CbNJSQkMHPmTMLCwoiMjCRPnjy8/PLLvPLKKzzwwAMAHPrl6utvHDwem5GlGuOVLLhNmjl58iSff/45H3/8Mbt376ZEiRIMHTqUjh07kj9//ku29c/vx4GrhLR/fr8rnjPGXMq6Sswt27VrF71796Z48eL07t2bYsWKMXXqVHbs2EGfPn2uCG2A1xqWxs83+yXP+flm57WGdoHSmJRYi9vcFFXl119/JSwsjB9//JHs2bPTpk0bevXqlaobQCVdgAyds42Dx2Pxz+/Haw1L24VJY1LBgtvckPPnz/Ptt98SFhbG2rVrKVSoEG+++SbdunXD39//ho7VvHIxC2pjboIFt0mVI0eOMGrUKEaMGMGff/5JuXLl+Oyzz2jbti1+ftYvbUxGsuA217V+/XqGDx/OpEmTiIuLo0mTJvTq1Yu6desiIk6XZ0yWZMFtrpCQkEBERARhYWEsXLiQ3Llz06lTJ3r27Enp0nbx0BinWXCbZKdOnWL8+PEMHz6cnTt3cvfdd/PBBx/QqVMnu3+IMR7EgtuwZ88ePvnkE8aOHcvJkycJDAzk/fffJyQkBB8f+0/EGE+Tqv9XikhRYKqq1hQRX+AHoCAwVlU/T88CTfpQVZYsWUJYWBg//PAD2bJlo1WrVvTq1YuqVas6XZ4x5jpSnIAjIgWAL4Gk5bRfAaJVNRBoKiL50rE+k8bi4uKYOHEiVapUoWbNmixcuJB+/fqxe/duvv76awttY7xAamZOJgBtgJPuv4OBKe7HS4ErZluISBcRiRaR6JiYmLSo09yimJgYBg0aRIkSJWjfvj1nz55l9OjR7Nu3j8GDB1O8eHGnSzTGpFKKXSWqehK4eOhXHuCA+/FJoOhV9hkDjAEICAjQtCjU3JyNGzcyfPhwJk6cyPnz52nUqBG9evWifv36NpzPGC91M1eeTgN+wAkgr/tv40ESExOZPXs2YWFhzJs3Dz8/P1588UV69uxJmTJlnC7PGHOLbia4VwNBwFSgIrA8TSsyN+306dN8+eWXDB8+nO3bt1O8eHGGDBlC586dKViwoNPlGWPSyM0E95fALBGpCZQFVqRtSeZG7d27l08++YTPPvuMEydOUK1aNb755htatGiBr6+v0+UZY9JYqoNbVYPd//uHiNTH1ep+S1UT0qm2LGX62gM3dKc8VWXZsmWEhYUxbdo0AFq2bEmvXr2oXr16RpVtjHHATc2uUNWD/DOyxNyi6WsP0H/aRmIvuL4DDxyPveb6i3FxcUydOpWwsDBWrVpFgQIF6Nu3L927d+fuu+/O8NqNMRnPpsV5gNA525JDO0nS+otJwf3XX38xZswYPv30Uw4ePEjp0qUZOXIk7du3J0+ePFc7rDEmk7Lg9gDXWmfx4PFYNm/eTFhYGBMmTODcuXM0aNCAsWPH0rBhQ7JlswWMjMmKLLg9wOXrL6omcm7XGuLW/0S5/0WTK1cuOnToQM+ePSlXrpyDlRpjPIEFtwd4rWFp+k/byKmjfxK7YwUnV88k/u/9FCxSlMGDB9O5c2cKFy7sdJnGGA9hwe2QxMREfvvtN5YsWUJUVBRHF0QSc2g/AHmKleaV9z5hSN8u5MiRw+FKjTGexoI7g5w7d45Vq1YRFRVFVFQUS5cu5fjx4wDcddddPF6zJjVq1KB27do8/PDDNh3dGHNNFtzp5OjRo8mt6SVLlhAdHU1cXBwAZcuWpXXr1gQFBREUFETJkiUtqI0xqWbBnQZUlV27diUHdVRUFFu2bAEgR44cVKlShd69exMUFMRjjz1GoUKFHK7YGOPNLLhvQnx8POvXr08O6aioKA4fPgxAgQIFCAwMpEOHDgQFBREQEECuXLkcrtgYk5lYcKfC6dOnWb58eXJIL1++nDNnzgBQsmRJ6tWrl9ztUaZMGRtfbYxJVxbcV3Ho0KHkvumoqCjWrVtHQkIC2bJlo2LFirz44osEBQVRo0YNW4DAGJPhsnxwqypbt269pNtj165dAPj5+VG9enXefPNNgoKCqF69OrfddpvDFRtjsrosF9znz59n9erVlwzLO3r0KAB33HEHQUFBdO/enaCgICpXrmy3RTXGeJxMH9zHjh1j2bJlyUG9cuVKzp8/D8CDDz5Is2bNkvunS5UqZcPyjDEeL1MFt6qyd+/eS7o9Nm3aBICPjw+PPvooPXr0ICgoiMDAQO644w6HKzbGmBvn1cGdkJDAxo0bk0N6yZIl7N/vmjaeL18+AgMDadOmDUFBQVStWpXcuXM7XLExxtw6rwrus2fPsnLlykv6p0+dOgVAsWLFqFmzZvJojwoVKpA9e3aHKzbGmLTnscE9fe0B3vt+OX9sXoNPzO/kPb6D3Vs3ER8fj4hQvnx52rZtm9w/fc8991j/tDEmS/DI4J6+9gD//uhr9k543fVEdl/8/EvTrH1XXmrZiMcee4wCBQo4W6QxxjjEI4M7dM42tFBJ8ge/QM5i5ch5ZynEx5dD+f1o3LiO0+UZY4yjPDK4Dx6PJVvO3NxereUVzxtjTFbnkTfV8M/vd0PPG2NMVuKRwf1aw9L4+V46IsTPNzuvNSztUEXGGOM5PLKrpHnlYoCrr/vg8Vj88/vxWsPSyc8bY0xW5pHBDa7wtqA2xpgr3XRXiYiME5GlIjIwLQsyxhhzfTcV3CLSAsiuqoGAv4g8kLZlGWOMuZabbXEHA1PcjxcAQRe/KCJdRCRaRKJjYmJuoTxjjDGXu9ngzgMccD8+CRS9+EVVHaOqAaoaUKRIkVupzxhjzGVuNrhPA0mDqvPewnGMMcbcoJsdVbIaV/fIcqAisO2aG65e/ZeI/HGT5wEoDPx1C/t7m6z2fsHec1Zh7/nGlLjWC6KqN3w0EbkNWAzMBxoB1VX1xE0Wl9K5olU1ID2O7Ymy2vsFe89Zhb3ntHNTXRyqehLXBcrlwOPpFdrGGGOudNMTcFT1GP+MLDHGGJNBvOGi4hinC8hgWe39gr3nrMLecxq5qT5uY4wxzvGGFrcxmYqIFBSR+iJS2OlajHey4PYwIlJURBY7XUdGEJHbRWS2iMwVkR9EJIfTNaU3EbkLiACqAgtFJMvMUHP/t73W6Toygoj4iMheEVnk/lchLY/vscGdFW9iJSIFgC9xzUzNCtoCw1S1PnAYeMLhejJCOaC3qr4HzAEecbiejDSUfybuZXYPA5NVNdj9b2NaHtwjgzsL38QqAWiD6zYCmZ6qjlDVue4/iwBHnKwnI6jqPFVdLiK1cLW6lzldU0YQkTrAGVxf0FlBdSBERKJEZJKIpOkttD0yuEnhJlaZlaqezIpj4kXkMaCAqi53upaMICKC6wv6Aq4v60zN3QX2FvCG07VkoFVAbVUNAo4DjdPy4J4a3Ne9iZXJPESkIPAJ8JLTtWQUdekOLAWaOl1PBngD+FRVjztdSAbaoKqH3I+3Amnaa+CpwW03scoC3C2xKUB/Vb2V+9l4DRHpJyId3H/mx9Uay+zqAd1FZBFQSUTGOlxPRpggIhVFJDsQAqxPy4N7aiAm3cQKXDex2uNcKSYddQQeBQa4r7y3cbqgDDAGaC8ivwLZgV8crifdqWqtpIt0wDpV7eR0TRngHWACsA5Ypqrz0vLgHjkBJyNvYmWMMd7GI4MbkofG1Qd+VdWsciXaGGNS5LHBbYwx5uo8tY/bGGPMNVhwG2OMl7HgNsYYL2PBbYwxXsaC2xhjvMz/A6bSIdGNm0JbAAAAAElFTkSuQmCC\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "import matplotlib.pyplot as plt\n",
    "from pylab import mpl\n",
    "import math\n",
    "\"\"\"一元多项式曲线拟合\n",
    "  拟合数据为xi=0, 1, 2, 3, 4, 5\n",
    "对应的函数值yi=2.1, 7.7, 13.6, 27.2, 40.9, 61\n",
    "\"\"\"\n",
    "x = [0, 1, 2, 3, 4, 5]\n",
    "y = [2.1, 7.7, 13.6, 27.2, 40.9, 61]\n",
    " \n",
    "\"\"\"完成拟合曲线参数计算前相应变量的计算\"\"\"\n",
    "def polynomial_fitting(data_x,data_y):\n",
    "    size=len(data_x)\n",
    "    i=0\n",
    "    sum_x = 0\n",
    "    sum_sqare_x =0\n",
    "    sum_third_power_x = 0\n",
    "    sum_four_power_x = 0\n",
    "    average_x = 0\n",
    "    average_y = 0\n",
    "    sum_y = 0\n",
    "    sum_xy = 0\n",
    "    sum_sqare_xy = 0\n",
    "    while i<size:\n",
    "        sum_x += data_x[i]\n",
    "        sum_y += data_y[i]\n",
    "        sum_sqare_x += math.pow(data_x[i],2)\n",
    "        sum_third_power_x +=math.pow(data_x[i],3)\n",
    "        sum_four_power_x +=math.pow(data_x[i],4)\n",
    "        sum_xy +=data_x[i]*data_y[i]\n",
    "        sum_sqare_xy +=math.pow(data_x[i],2)*data_y[i]\n",
    "        i += 1;\n",
    "    average_x=sum_x/size\n",
    "    average_y=sum_y/size\n",
    "    return [[size, sum_x, sum_sqare_x, sum_y]\n",
    "        , [sum_x, sum_sqare_x, sum_third_power_x, sum_xy]\n",
    "        , [sum_sqare_x,sum_third_power_x,sum_four_power_x,sum_sqare_xy]]\n",
    " \n",
    " \n",
    "\"\"\"完成拟合曲线参数的计算\n",
    " 其中解方程的时候，利用高斯消元法计算相应的参数值\n",
    "\"\"\"\n",
    "def calculate_parameter(data):\n",
    "    #i用来控制列元素，line是一行元素,j用来控制循环次数,datas用来存储局部变量。保存修改后的值\n",
    "    i = 0;\n",
    "    j = 0;\n",
    "    line_size = len(data)\n",
    " \n",
    "   #将行列式变换为上三角行列式\n",
    "    while j < line_size-1:\n",
    "        line = data[j]\n",
    "        temp = line[j]\n",
    "        templete=[]\n",
    "        for x in line:\n",
    "            x=x/temp\n",
    "            templete.append(x)\n",
    "        data[j]=templete\n",
    "        #flag标志应该进行消元的行数\n",
    "        flag = j+1\n",
    "        while flag < line_size:\n",
    "            templete1 = []\n",
    "            temp1=data[flag][j]\n",
    "            i = 0\n",
    "            for x1 in data[flag]:\n",
    "                if x1!=0:\n",
    "                   x1 = x1-(temp1*templete[i])\n",
    "                   templete1.append(x1)\n",
    "                else:\n",
    "                   templete1.append(0)\n",
    "                i += 1\n",
    "            data[flag] = templete1\n",
    "            flag +=1\n",
    "        j += 1\n",
    " \n",
    " \n",
    "    #求相应的参数值\n",
    " \n",
    "    parameters=[]\n",
    "    i=line_size-1\n",
    "    #j标识减去元素个数\n",
    "    #flag_rol标识除那一列\n",
    "    flag_j=0\n",
    "    rol_size=len(data[0])\n",
    "    flag_rol=rol_size-2\n",
    "    #获得解的个数\n",
    "    while i>=0:\n",
    "        operate_line = data[i]\n",
    "        if i==line_size-1:\n",
    "            parameter=operate_line[rol_size-1]/operate_line[flag_rol]\n",
    "            parameters.append(parameter)\n",
    "        else:\n",
    "            flag_j=(rol_size-flag_rol-2)\n",
    "            temp2=operate_line[rol_size-1]\n",
    "            #result_flag为访问已求出解的标志\n",
    "            result_flag=0\n",
    "            while flag_j>0:\n",
    "                temp2-=operate_line[flag_rol+flag_j]*parameters[result_flag]\n",
    "                result_flag+=1\n",
    "                flag_j-=1\n",
    "            parameter=temp2/operate_line[flag_rol]\n",
    "            parameters.append(parameter)\n",
    "        flag_rol-=1\n",
    "        i-=1\n",
    "    return parameters\n",
    " \n",
    "\"\"\"计算拟合曲线的值\"\"\"\n",
    "def calculate(data_x,parameters):\n",
    "    datay=[]\n",
    "    for x in data_x:\n",
    "        datay.append(parameters[2]+parameters[1]*x+parameters[0]*x*x)\n",
    "    return datay\n",
    " \n",
    " \n",
    " \n",
    "\"\"\"完成函数的绘制\"\"\"\n",
    "def draw(data_x,data_y_new,data_y_old):\n",
    "    plt.plot(data_x,data_y_new,label=\"拟合曲线\",color=\"black\")\n",
    "    plt.scatter(data_x,data_y_old,label=\"离散数据\")\n",
    "    mpl.rcParams['font.sans-serif'] = ['SimHei']\n",
    "    mpl.rcParams['axes.unicode_minus'] = False\n",
    "    plt.title(\"一元多项式拟合数据\")\n",
    "    plt.legend(loc=\"upper left\")\n",
    "    plt.show()\n",
    " \n",
    "data=polynomial_fitting(x,y)\n",
    "parameters=calculate_parameter(data)\n",
    "for w in parameters:\n",
    "    print(w)\n",
    "newData=calculate(x,parameters)\n",
    "draw(x,newData,y)\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "id": "45591581",
   "metadata": {},
   "outputs": [
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "C:\\Users\\Administrator\\AppData\\Local\\Temp\\ipykernel_1612\\3156658090.py:120: MatplotlibDeprecationWarning: Axes3D(fig) adding itself to the figure is deprecated since 3.4. Pass the keyword argument auto_add_to_figure=False and use fig.add_axes(ax) to suppress this warning. The default value of auto_add_to_figure will change to False in mpl3.5 and True values will no longer work in 3.6.  This is consistent with other Axes classes.\n",
      "  ax = Axes3D(fig)\n"
     ]
    },
    {
     "ename": "ValueError",
     "evalue": "shape mismatch: objects cannot be broadcast to a single shape.  Mismatch is between arg 0 with shape (5,) and arg 1 with shape (6,).",
     "output_type": "error",
     "traceback": [
      "\u001b[1;31m---------------------------------------------------------------------------\u001b[0m",
      "\u001b[1;31mValueError\u001b[0m                                Traceback (most recent call last)",
      "Input \u001b[1;32mIn [2]\u001b[0m, in \u001b[0;36m<cell line: 132>\u001b[1;34m()\u001b[0m\n\u001b[0;32m    130\u001b[0m data2 \u001b[38;5;241m=\u001b[39m calculate_parameter(data)\n\u001b[0;32m    131\u001b[0m newY \u001b[38;5;241m=\u001b[39m calculate(x1, x2, data2)\n\u001b[1;32m--> 132\u001b[0m \u001b[43mdraw\u001b[49m\u001b[43m(\u001b[49m\u001b[43m[\u001b[49m\u001b[43mx1\u001b[49m\u001b[43m,\u001b[49m\u001b[43m \u001b[49m\u001b[43mx2\u001b[49m\u001b[43m]\u001b[49m\u001b[43m,\u001b[49m\u001b[43m \u001b[49m\u001b[43my\u001b[49m\u001b[43m,\u001b[49m\u001b[43m \u001b[49m\u001b[43mnewY\u001b[49m\u001b[43m)\u001b[49m\n",
      "Input \u001b[1;32mIn [2]\u001b[0m, in \u001b[0;36mdraw\u001b[1;34m(data_x, old_y, new_y)\u001b[0m\n\u001b[0;32m    119\u001b[0m \u001b[38;5;66;03m#创建Axes3D对象，让其包含图像3D坐标\u001b[39;00m\n\u001b[0;32m    120\u001b[0m ax \u001b[38;5;241m=\u001b[39m Axes3D(fig)\n\u001b[1;32m--> 121\u001b[0m \u001b[43max\u001b[49m\u001b[38;5;241;43m.\u001b[39;49m\u001b[43mscatter\u001b[49m\u001b[43m(\u001b[49m\u001b[43mdata_x\u001b[49m\u001b[43m[\u001b[49m\u001b[38;5;241;43m0\u001b[39;49m\u001b[43m]\u001b[49m\u001b[43m,\u001b[49m\u001b[43m \u001b[49m\u001b[43mdata_x\u001b[49m\u001b[43m[\u001b[49m\u001b[38;5;241;43m1\u001b[39;49m\u001b[43m]\u001b[49m\u001b[43m,\u001b[49m\u001b[43m \u001b[49m\u001b[43mold_y\u001b[49m\u001b[43m,\u001b[49m\u001b[43m \u001b[49m\u001b[43mcolor\u001b[49m\u001b[38;5;241;43m=\u001b[39;49m\u001b[38;5;124;43m'\u001b[39;49m\u001b[38;5;124;43mred\u001b[39;49m\u001b[38;5;124;43m'\u001b[39;49m\u001b[43m)\u001b[49m\n\u001b[0;32m    122\u001b[0m ax\u001b[38;5;241m.\u001b[39mplot(data_x[\u001b[38;5;241m0\u001b[39m], data_x[\u001b[38;5;241m1\u001b[39m], new_y, color\u001b[38;5;241m=\u001b[39m\u001b[38;5;124m'\u001b[39m\u001b[38;5;124mblack\u001b[39m\u001b[38;5;124m'\u001b[39m)\n\u001b[0;32m    123\u001b[0m mpl\u001b[38;5;241m.\u001b[39mrcParams[\u001b[38;5;124m'\u001b[39m\u001b[38;5;124mfont.sans-serif\u001b[39m\u001b[38;5;124m'\u001b[39m] \u001b[38;5;241m=\u001b[39m [\u001b[38;5;124m'\u001b[39m\u001b[38;5;124mSimHei\u001b[39m\u001b[38;5;124m'\u001b[39m]\n",
      "File \u001b[1;32mc:\\users\\administrator\\appdata\\local\\programs\\python\\python38\\lib\\site-packages\\matplotlib\\__init__.py:1412\u001b[0m, in \u001b[0;36m_preprocess_data.<locals>.inner\u001b[1;34m(ax, data, *args, **kwargs)\u001b[0m\n\u001b[0;32m   1409\u001b[0m \u001b[38;5;129m@functools\u001b[39m\u001b[38;5;241m.\u001b[39mwraps(func)\n\u001b[0;32m   1410\u001b[0m \u001b[38;5;28;01mdef\u001b[39;00m \u001b[38;5;21minner\u001b[39m(ax, \u001b[38;5;241m*\u001b[39margs, data\u001b[38;5;241m=\u001b[39m\u001b[38;5;28;01mNone\u001b[39;00m, \u001b[38;5;241m*\u001b[39m\u001b[38;5;241m*\u001b[39mkwargs):\n\u001b[0;32m   1411\u001b[0m     \u001b[38;5;28;01mif\u001b[39;00m data \u001b[38;5;129;01mis\u001b[39;00m \u001b[38;5;28;01mNone\u001b[39;00m:\n\u001b[1;32m-> 1412\u001b[0m         \u001b[38;5;28;01mreturn\u001b[39;00m \u001b[43mfunc\u001b[49m\u001b[43m(\u001b[49m\u001b[43max\u001b[49m\u001b[43m,\u001b[49m\u001b[43m \u001b[49m\u001b[38;5;241;43m*\u001b[39;49m\u001b[38;5;28;43mmap\u001b[39;49m\u001b[43m(\u001b[49m\u001b[43msanitize_sequence\u001b[49m\u001b[43m,\u001b[49m\u001b[43m \u001b[49m\u001b[43margs\u001b[49m\u001b[43m)\u001b[49m\u001b[43m,\u001b[49m\u001b[43m \u001b[49m\u001b[38;5;241;43m*\u001b[39;49m\u001b[38;5;241;43m*\u001b[39;49m\u001b[43mkwargs\u001b[49m\u001b[43m)\u001b[49m\n\u001b[0;32m   1414\u001b[0m     bound \u001b[38;5;241m=\u001b[39m new_sig\u001b[38;5;241m.\u001b[39mbind(ax, \u001b[38;5;241m*\u001b[39margs, \u001b[38;5;241m*\u001b[39m\u001b[38;5;241m*\u001b[39mkwargs)\n\u001b[0;32m   1415\u001b[0m     auto_label \u001b[38;5;241m=\u001b[39m (bound\u001b[38;5;241m.\u001b[39marguments\u001b[38;5;241m.\u001b[39mget(label_namer)\n\u001b[0;32m   1416\u001b[0m                   \u001b[38;5;129;01mor\u001b[39;00m bound\u001b[38;5;241m.\u001b[39mkwargs\u001b[38;5;241m.\u001b[39mget(label_namer))\n",
      "File \u001b[1;32mc:\\users\\administrator\\appdata\\local\\programs\\python\\python38\\lib\\site-packages\\mpl_toolkits\\mplot3d\\axes3d.py:2379\u001b[0m, in \u001b[0;36mAxes3D.scatter\u001b[1;34m(self, xs, ys, zs, zdir, s, c, depthshade, *args, **kwargs)\u001b[0m\n\u001b[0;32m   2376\u001b[0m had_data \u001b[38;5;241m=\u001b[39m \u001b[38;5;28mself\u001b[39m\u001b[38;5;241m.\u001b[39mhas_data()\n\u001b[0;32m   2377\u001b[0m zs_orig \u001b[38;5;241m=\u001b[39m zs\n\u001b[1;32m-> 2379\u001b[0m xs, ys, zs \u001b[38;5;241m=\u001b[39m \u001b[43mnp\u001b[49m\u001b[38;5;241;43m.\u001b[39;49m\u001b[43mbroadcast_arrays\u001b[49m\u001b[43m(\u001b[49m\n\u001b[0;32m   2380\u001b[0m \u001b[43m    \u001b[49m\u001b[38;5;241;43m*\u001b[39;49m\u001b[43m[\u001b[49m\u001b[43mnp\u001b[49m\u001b[38;5;241;43m.\u001b[39;49m\u001b[43mravel\u001b[49m\u001b[43m(\u001b[49m\u001b[43mnp\u001b[49m\u001b[38;5;241;43m.\u001b[39;49m\u001b[43mma\u001b[49m\u001b[38;5;241;43m.\u001b[39;49m\u001b[43mfilled\u001b[49m\u001b[43m(\u001b[49m\u001b[43mt\u001b[49m\u001b[43m,\u001b[49m\u001b[43m \u001b[49m\u001b[43mnp\u001b[49m\u001b[38;5;241;43m.\u001b[39;49m\u001b[43mnan\u001b[49m\u001b[43m)\u001b[49m\u001b[43m)\u001b[49m\u001b[43m \u001b[49m\u001b[38;5;28;43;01mfor\u001b[39;49;00m\u001b[43m \u001b[49m\u001b[43mt\u001b[49m\u001b[43m \u001b[49m\u001b[38;5;129;43;01min\u001b[39;49;00m\u001b[43m \u001b[49m\u001b[43m[\u001b[49m\u001b[43mxs\u001b[49m\u001b[43m,\u001b[49m\u001b[43m \u001b[49m\u001b[43mys\u001b[49m\u001b[43m,\u001b[49m\u001b[43m \u001b[49m\u001b[43mzs\u001b[49m\u001b[43m]\u001b[49m\u001b[43m]\u001b[49m\u001b[43m)\u001b[49m\n\u001b[0;32m   2381\u001b[0m s \u001b[38;5;241m=\u001b[39m np\u001b[38;5;241m.\u001b[39mma\u001b[38;5;241m.\u001b[39mravel(s)  \u001b[38;5;66;03m# This doesn't have to match x, y in size.\u001b[39;00m\n\u001b[0;32m   2383\u001b[0m xs, ys, zs, s, c \u001b[38;5;241m=\u001b[39m cbook\u001b[38;5;241m.\u001b[39mdelete_masked_points(xs, ys, zs, s, c)\n",
      "File \u001b[1;32m<__array_function__ internals>:180\u001b[0m, in \u001b[0;36mbroadcast_arrays\u001b[1;34m(*args, **kwargs)\u001b[0m\n",
      "File \u001b[1;32mc:\\users\\administrator\\appdata\\local\\programs\\python\\python38\\lib\\site-packages\\numpy\\lib\\stride_tricks.py:539\u001b[0m, in \u001b[0;36mbroadcast_arrays\u001b[1;34m(subok, *args)\u001b[0m\n\u001b[0;32m    532\u001b[0m \u001b[38;5;66;03m# nditer is not used here to avoid the limit of 32 arrays.\u001b[39;00m\n\u001b[0;32m    533\u001b[0m \u001b[38;5;66;03m# Otherwise, something like the following one-liner would suffice:\u001b[39;00m\n\u001b[0;32m    534\u001b[0m \u001b[38;5;66;03m# return np.nditer(args, flags=['multi_index', 'zerosize_ok'],\u001b[39;00m\n\u001b[0;32m    535\u001b[0m \u001b[38;5;66;03m#                  order='C').itviews\u001b[39;00m\n\u001b[0;32m    537\u001b[0m args \u001b[38;5;241m=\u001b[39m [np\u001b[38;5;241m.\u001b[39marray(_m, copy\u001b[38;5;241m=\u001b[39m\u001b[38;5;28;01mFalse\u001b[39;00m, subok\u001b[38;5;241m=\u001b[39msubok) \u001b[38;5;28;01mfor\u001b[39;00m _m \u001b[38;5;129;01min\u001b[39;00m args]\n\u001b[1;32m--> 539\u001b[0m shape \u001b[38;5;241m=\u001b[39m \u001b[43m_broadcast_shape\u001b[49m\u001b[43m(\u001b[49m\u001b[38;5;241;43m*\u001b[39;49m\u001b[43margs\u001b[49m\u001b[43m)\u001b[49m\n\u001b[0;32m    541\u001b[0m \u001b[38;5;28;01mif\u001b[39;00m \u001b[38;5;28mall\u001b[39m(array\u001b[38;5;241m.\u001b[39mshape \u001b[38;5;241m==\u001b[39m shape \u001b[38;5;28;01mfor\u001b[39;00m array \u001b[38;5;129;01min\u001b[39;00m args):\n\u001b[0;32m    542\u001b[0m     \u001b[38;5;66;03m# Common case where nothing needs to be broadcasted.\u001b[39;00m\n\u001b[0;32m    543\u001b[0m     \u001b[38;5;28;01mreturn\u001b[39;00m args\n",
      "File \u001b[1;32mc:\\users\\administrator\\appdata\\local\\programs\\python\\python38\\lib\\site-packages\\numpy\\lib\\stride_tricks.py:421\u001b[0m, in \u001b[0;36m_broadcast_shape\u001b[1;34m(*args)\u001b[0m\n\u001b[0;32m    416\u001b[0m \u001b[38;5;124;03m\"\"\"Returns the shape of the arrays that would result from broadcasting the\u001b[39;00m\n\u001b[0;32m    417\u001b[0m \u001b[38;5;124;03msupplied arrays against each other.\u001b[39;00m\n\u001b[0;32m    418\u001b[0m \u001b[38;5;124;03m\"\"\"\u001b[39;00m\n\u001b[0;32m    419\u001b[0m \u001b[38;5;66;03m# use the old-iterator because np.nditer does not handle size 0 arrays\u001b[39;00m\n\u001b[0;32m    420\u001b[0m \u001b[38;5;66;03m# consistently\u001b[39;00m\n\u001b[1;32m--> 421\u001b[0m b \u001b[38;5;241m=\u001b[39m \u001b[43mnp\u001b[49m\u001b[38;5;241;43m.\u001b[39;49m\u001b[43mbroadcast\u001b[49m\u001b[43m(\u001b[49m\u001b[38;5;241;43m*\u001b[39;49m\u001b[43margs\u001b[49m\u001b[43m[\u001b[49m\u001b[43m:\u001b[49m\u001b[38;5;241;43m32\u001b[39;49m\u001b[43m]\u001b[49m\u001b[43m)\u001b[49m\n\u001b[0;32m    422\u001b[0m \u001b[38;5;66;03m# unfortunately, it cannot handle 32 or more arguments directly\u001b[39;00m\n\u001b[0;32m    423\u001b[0m \u001b[38;5;28;01mfor\u001b[39;00m pos \u001b[38;5;129;01min\u001b[39;00m \u001b[38;5;28mrange\u001b[39m(\u001b[38;5;241m32\u001b[39m, \u001b[38;5;28mlen\u001b[39m(args), \u001b[38;5;241m31\u001b[39m):\n\u001b[0;32m    424\u001b[0m     \u001b[38;5;66;03m# ironically, np.broadcast does not properly handle np.broadcast\u001b[39;00m\n\u001b[0;32m    425\u001b[0m     \u001b[38;5;66;03m# objects (it treats them as scalars)\u001b[39;00m\n\u001b[0;32m    426\u001b[0m     \u001b[38;5;66;03m# use broadcasting to avoid allocating the full array\u001b[39;00m\n",
      "\u001b[1;31mValueError\u001b[0m: shape mismatch: objects cannot be broadcast to a single shape.  Mismatch is between arg 0 with shape (5,) and arg 1 with shape (6,)."
     ]
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAATwAAAE1CAYAAAB6ClzOAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjUuMSwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/YYfK9AAAACXBIWXMAAAsTAAALEwEAmpwYAAB6cklEQVR4nO29eXxb5Zk2fB1Jlmx53+N4d5zFseMlsRMoLYWBoSUsCW9ZQjssTVt4p+20fHS6TVo++rXTdtoynXnpQqdAaZm30JLAJEBIp7QDtGxJIPG+7/Km3ZIlWevz/WGeh6Nj2dZyJB3b5/r9+AG2fM6xrHOd+3nu674ujhACGTJkyNgMUCT7AmTIkCEjUZAJT4YMGZsGMuHJkCFj00AmPBkyZGwayIQnQ4aMTQPVGt+XW7gyZMhIFjixDyhXeDJkyNg0kAlPhgwZmwYy4cmQIWPTQCY8GTJkbBrIhCdDhoxNA5nwZMiQsWkgE54MGTI2DWTCkyFDxqaBTHgyZMjYNJAJT4YMGZsGMuHJkCFj00AmPBkyZGwayIQnQ4aMTQOZ8GTIkLFpIBOeDBkyNg1kwpMhQ8amgUx4MmTI2DSQCU+GDBmbBjLhyZAhY9NAJjwZMmRsGsiEJ0OGjE2DtVLLZGwAEELgdrsRCASQkpICpVIJhUJ+1snYfOAIWTWJUY5pXOcIBALweDzwer3w+XzgOA4GgwG5ublIT0+HSqWSCVCGVCF6TKNc4W1QEELg8/kwODiIzMxMaDQaAEBaWhrMZjO0Wi1UKhXcbjcAQKlUsupPpVKB40T/rMmQkXTIhLcBQQiBx+NBIBBAIBDAzMwMnE4nlEolPB4P/H4/tFotNBoNNBoNCCEIBAJYXFwEIQQcxzECpBWgTIAyNgLkJe0Gg8/ng9frBQB4PB6cPXsW6enp2LNnDwghIISgs7MTqampcDqd8Pl8yMrKQm5uLnJycqBWq9nr+J8NmQBlJAHyklZGaNAlLN2nM5lM6O/vR15eHvLy8qBQKOD3+6FQKKDRaFBSUoKsrCwEAgHMz8/DarVCp9PB7/cjOzubEWBKSkpQBUghE6CM9QiZ8DYAAoEAvF4vAoEACCEYHByE3W5Ha2srpqenl72e4zhWvSkUCuTm5iI3NxfV1dXw+/2Yn5+HxWLBxMQECCFBBKhSqRgBulwuRnQyAcpYD5AJbx2DEAK/38+WsIuLi+js7ERhYSH27dsHjuOCyI3uz60GpVLJqkJgaYlMCXB8fByEEOTk5CA3NxfZ2dlQKpUhCVClUrF/ZAKUIRXIhLdOQQiB1+uF3+8Hx3GYnZ3FyMgI6uvrkZOTw17HJ7zVvrYSVCoV8vPzkZ+fD2CJAK1WK8xmM0ZHR8FxXEgCpETMrwCVSiXUarVMgDKSBpnw1iGoto5WVn19ffD5fNi/fz9SUlKCXhsr4QmhUqlQUFCAgoICAIDX64XVaoXRaMTw8DCUSmUQASoUChBCYLVaMTMzgx07drDj0CWwQqGQCVBGQiAT3jqCsDHhcDjQ2dmJiooKlJaWhiSNlcgtWsITIiUlBYWFhSgsLASw1Bm2Wq3Q6/UYGhqCSqVCbm4u1Go1AARVgD6fj10jfwksE6CMeEEmvHUCvrYOACYnJzE1NYXGxkZkZGSs+HMcx7Gf4X8tXlCr1SgqKkJRUREAwO12w2q1Ym5uDlarFW63mzVJMjIyWAXIl9PIBCgjXpAJbx2A7odRYuju7oZarcb+/fuhVCpX/Vmxl7SRQqPRoLi4GGlpaZienkZVVRUsFgumpqZgt9uh0WiCCJBem9frXUaA/DlgmQBlRAOZ8CQM4RJ2fn4ePT092LZtG7Zs2RLWMeK9pI0UqampKCkpQUlJCQDA5XLBYrFgcnISCwsLSE1NZQSYnp6+IgHyJTAyAcoIFzLhSRR8bR0AjIyMwGg0Yu/evUhLS4v4eHxCkBI5pKWlIS0tDVu3bgUhhBHg2NgYHA4HtFotI0CtVssI0OPxwO12M+kNJUA6Byyl31GGdCATnsRAN/R1Oh3cbjeKi4vR1dWFnJwctLW1RexqkuwlbSTgOA5arRZarRalpaUghMDpdMJisWBkZAROpxPp6emMANPS0lYlQLoElglQBoVMeBICX1sHADabDVNTU9i1axfTwUUKqS1pIwHHcUhPT0d6ejrKyspACIHD4YDFYsHQ0BAWFxeRkZHBZDBpaWns9/J4POwfp9OJoqIi2QtQhkx4UgFfW0cIwfT0NBwOBy655BIm6YgGK1V46xEcxyEjIwMZGRkoLy8HIQQLCwuwWCwYGBiA2+1GZmYmI8DU1FS43W4YjUbk5OTA4/EAWBqnEzZBZGwOyISXZAinElwuFzo7O5Geno7MzMyYyA5YX0vaSMFxHDIzM5GZmYmKigoEAgFGgH19ffB4PEhNTUUgEIDP52NWWMD7FSAgE+Bmgkx4SQRfW8dxHGZmZjA2Nob6+np4vV6YTKaYz7GRCU8IhUKBrKwsZGVlobKyEoFAANPT05idnUVPTw+8Xm+QEQLfCktIgMIusIyNAZnwkgT+Etbv96O3txcAsH//fqhUKphMJtFIab3u4cUKhUKB9PR0ZGdnY/v27cwKy2KxrGqFRTNAZDfojQeZ8BIMobZuYWEBXV1dqKysRGlpKXsdnUCIFaGqk810w/IdYvhWWAAissKSvQA3BmTCSyCE2rrx8XHMzMygqakJ6enpQa8NNRIWy3mFx94MFd5aWM0Ka2xsDACCjBBkAlz/kAkvARD61nm9XnR3dyMtLQ0HDhxYsQoTg5RWuvk2C+GF4wFIEa0VlkyA6wcy4cUZwiWsxWJBb28vtm/fzgbsQ0GhUIhS4VHilOqkRbwRCeEJEa0VFjVDnZqaQkFBAdLT02UClAhkwosjhNbrw8PDsFgs2LdvH1JTU1f9WTErPEIIy63Iy8sTdbksdcRCeEKEa4WVm5uLzMxMLCwsoKCgQLbDlxBkwosDhNo6t9uNzs5O5OXloa2tLawPuJj7bDabDXa7HXl5eRgYGMDCwgI0Gg1SUlKCvOpkRIZQVlgWiwUzMzPo7++H2+2GRqNBUVFRkBVWKDt82gWWCTC+kAlPZAi1dfTpX1dXxzbHw4EYhOfxeNDf3w+/349LLrkEPp8PCoUC09PTWFhYgNPpxNTUFPx+P1ua0e7kRoGYFd5a0Gg02LJlC3OyuXDhAjQazZpWWHwzVEB2g44nNs4nWwIIBAIYGxtDSUkJAoEABgYGsLi4iLa2toirqFj38KxWK7q7u1FaWgqbzRbUGFEoFEhNTUVFRQVLKrNaraw7yXEcuzHp3tR6RSIJTwiO41BcXIzy8nIA71thTUxMYGFhAWlpaSGtsGQ36PhBJjwRwG9MjI+PIzs7G52dndi6dSvq6uqi+oBGW+ERQjA+Po7Z2Vm0tLTA7/fDZrOFfB2FUqkM6k56vV5YLBbo9XoMDg5CrVYH7U2tpxsumYQnPHe0VliyG7R4kAkvRgit171eL9rb27Fnzx5kZWVFfdxoCM/r9aKrqwupqanYv38/FAoFFhYWIjYPSElJCdqbWlxcZNMJdrs95I0pYzkCgcCK1XEsVlgyAUYPmfBiAP+D5/f70dPTA7/fjwMHDqzZhV0LkS5p5+fn0dXVtcwNWYxZWr5LcagbMzMzE7m5ucjLy4NGown7uImAlCq81bCWFZbL5UJGRkYQAdJzyHb44UMmvCgg1NbZbDZ0d3ejqqoKXq9XlD2vcEmJEIKJiQnMzMygpaUFWq12zePE0hAJdWPa7XZYLJZlA/q5ubnLYiMTjfVCeEJEY4VFz+n1eqHT6ZCSkoL8/HzZDZoHmfAihHA8bGxsDHNzc2huboZWq8Xs7KwoGrdwZmnpElaj0bAlrBB8covHB53juGUOJXQ8a3JyEoQQ5OTkQKPRbBrtHx9ivefhWGFlZWUxAnS73UhJSQHHccwNGnjfCWazukHLhBcmQo2HdXV1IT09PYhsxJqQoOdcCfPz8+ju7kZ1dTULxAmFRNtDCQf06XjW7OwsLBYL3n33Xbb8zczMjHsHOJkVXjwRygrLZrMxHaDD4UBmZiZ74FCVQCgrrM3kBSgTXhjgW69zHAez2Yy+vj7s2LGDqe4pxJpiWG0GlmbShjIdCHWcZPrh0fEstVoNlUqFmpoaWCwWTE9PB2nT8vLymDRDTBBCNvxNDCwRV05ODnJycgAAAwMDSE1NxcLCwqpWWMDmMkOVCW8NCK3Xh4aGMD8/v+J4mJgVnhA+nw9dXV1ISUkJK5MWkJ4zilqtRnFxMYqLiwG8r00bHx/HwsJCyM6kjMhBra6ys7MBhG+FBWxsApQJbwUIGxOLi4vo7OxEQUEBWltbV7wR40V4NpsNXV1dqKqqwtatW2M6lpRIUKhNE4b00A5wbm5uVB3gjbqkXQtCSUy0VlgbzQ1aJrwQEI6Hzc3NYXh4GLt372Z7UytBLONO/rXodDrodDo0NjYiIyMjop9P9pI2Egg7k4FAIKgD7PP5gjrA4YzAyYQXGtFaYQndoPkEuB7coGXCE4A2JuiQd39/PzweD/bv3x+WxELspkVHRweUSmXYS1ghQtlDrRcoFAq2LKuqqgpalo2Pj4e8KWUsIRAIRPR+RGuFtZIdvlSdYGTCew/CJazD4UBXVxfKyspQVlYWkYBUDMKz2+1wOByoqqoKsn6PFOupwlsLwmWZ8Kak9ky0A7yeyT5W+P3+mJabkVph8Z1gFhcXMT4+ji1btkCr1cJqtcLpdGLnzp1i/XpRQyY8LNfW0SVkQ0MDMjMzIzpWrBUeIQRTU1OYnJxEWlpaTGQHbCzCE0J4U1J7JjoCR4O5c3NzE058yX5/11rSRoq1rLD489YZGRlYWFhgD/9XX30VfX19+M53viPa9USLTU14VFs3MjKC8vJy+Hw+9PT0RNQFFSIWwqPn5zgObW1tOHv2bFTH4WMlwtuI4Nsz0eH8wcFB6PV6zMzMsNGsvLy8mEf/1kKyK0uxCU8IoRUWf956YWEBbrcbs7Oz8Pv9WFhYiHjv+ejRo/jVr36lB6AnhDQIv88tvbn/DuAgACeAuwkh76513E1LeHxt3dTUFLKzs9Hd3b1sFjVSREt4drsdnZ2dqKioQFlZWdTnF4J/0/FvwmRXIPEGHc7PyspCRkYGCgoK2GQCNefMyspCXl5ekDBXLCSb8BKtP+TPWwPA2bNnoVar8cMf/hAvv/wycnJyoNVqceWVV6K+vn7Na7v77rvxq1/96qMAfrPCS64FsP29fw4A+Pl7/14Vm5Lw+No6YKk87+/vDzmLGimi6dJOTU1hfHwce/bsWbaEFvPGocfZKEvacEDfv1CjWXQyQafTIRAIIDs7G3l5eUyWEQviXWGtB5SWluLf/u3f8G//9m9ISUlBamoqvve97+GRRx5Z00no8ssvBwDzKi85BOA3ZOmD/BbHcTkcx5UQQmZWO+6mIjyh9brX60VHRwcIIWhraxPlAxpJhUcdVgghLICbj3htum9GwhOCP5lQXV0dpEsbHR1l38/Ly0NWVlbEn41kV3hSgsvlwq5du3DLLbfgnnvuEeuwpQAmef+ve+9rMuEBy7V1RqMRAwMD2LlzJwYHB0U7T7hd2oWFBXR2dqK8vBylpaUhb47NREzxRDjEI9SlURPUubk5DAwMsE35vLw8Zs++GjbLSNtK4L8/CwsLETf/4oVNQXjC8bDBwUHY7Xa0trZCo9FgeHhYtCVIOBXe9PQ0xsbGQi5h+RCb8Kg2azMRabS/50omqNSenZqg5uXlhRyBow9WGYDT6Yy4aREGpgCU8/6/7L2vrYoNTXihxsM6OjpQXFyMffv2BcXmiSUWXo3w/H4/ent74fP5Qi5hIzlWJCCEoKenB0ajERqNBhkZGaza3ehViFhLy5VMUKk5p9AEdTNXeMKHjMPhWNPkIgqcAvB5juOexlKzYn6t/TtgAxMeX1vHcRxmZ2cxOjqK+vp6NlBNoVAo4Pf7RTnvSiTlcDjQ0dGB0tJSlJeXJyyq0eVywel0YuvWrThw4AB8Ph9mZ2dhNBpx/vx5pKamMjHvRhzWj9ceaCgTVLPZzEbg0tLS4PF44PV6E26CmmzfQb/fHyTpolZVkeD2228HgDcBFHAcpwPw/wJIAQBCyCMATmNJkjKEJVnKJ8M57oYjPKFvnd/vR19fHwKBwIpVlZjjYKG6tDMzMxgdHUVDQ0NEORexEp7RaER/fz80Gg2qqqrg8XigVqtRWFiI+fl5NDQ0wOVywWw2Bw3r5+Xlbai82niTON8ElY7Azc7OYmpqijXFaAMkESNwhJCkjtmFIrxIl7RPPfUUnnrqqRWNHt/rzn4u0mvbUIQnXMLSxkBFRcWKjQEgfktaSrZerzesJexqx4oEhBCMjIzAbDajtbUV77zzTtD3+d1fGiRTVlbGhvXNZjOmpqYQCAQSeqPGA8nYq1QqlcjIyEB2djZ27tzJBvP5c6nxNEGNdaxMjPPHWuHFCxuG8ITjYRMTE5ieng7LYUTMJS3t0tIl7NatW1FRUZGwqEYqtcnIyMC+ffsi+uDzh/WpVEM4q0qXv+F0KqWAZMlD+OcVDuZ7PB6YzeYgE1RaVYthgprsvVkh4blcrrhPtoSLdU94Qm2dz+dDd3c3y3kIpyoRe0nrcDjQ3t4ecr8wEkRKeDabDZ2dnaitrWUGm7EcV3ijCjuV6enpjACl8oEWIlmEtxrpqNXqoLEsuq1A82np+8pPJxPr3ImAkPAASKaBs64JT6its1qt6OnpWfOGF0KsJW0gEMDo6CicTicuu+yymDerI5na0Ol0mJycRHNzc1S27+FA2Kl0OBzM7t7j8bBJhXC96jYyIiFaahJB82mpCSo/nSySfVUpEB49v9QE2Ov2UykcDxseHobZbMbevXsjfiqKUeE5nU50dHQgPz8fgUBAlM5cOCJmKnWhTZlE2b7zzTrpqNb8/DzMZjPzqqP7VLEEkseKZC5poyGdlUxQ6b6q3+9nvnTUml2ISL3wxIawwpMS6a07wqMuGMASUXk8HnR2diInJwetra1Rfchi3cObm5vD0NAQ6uvrkZqaip6enqiPxcdaxERJduvWrWFLXeIFYVoZnVSYnZ3FwMAAc8x1OBzQarUJu9ZkLmnFOK9wX5WaoNIlMDVBpY0l+llOdoVHCU9qAux1RXh0CdvR0YHt27fD6XRicHAQu3btYiNB0SDaCo86IrtcLuaITJfYYmC16zIYDBgYGEB9fT1LqgoXiZi0EE4q6PV6TE5OYmRkZJlQd6PIX/iI17JyJRNUvjEnzf5IFtnzCc/lcsVsyCEm1g3h+Xw+pq1TKBQYHh6Gz+dDW1tbzDeMUqmMuMJzOp3o7OxEUVERdu3aFeREIhbhhSImfnJatL97MkbLUlNTkZ6ejl27doEQwpxKurq64Pf7WXWYk5Mj6nJMCl3aeCKUCer4+DisVivOnj2LtLQ0tv+XqMo6EAiwz2WcpiyihuQJT6itczqdMJlMKCkpQWNjo2jLBkqm4UCv12NwcDBkdSVmx1dITLS6zc7ODhqNCwdSmp3lOC4oq4IfIDMyMhLSqj1abHTCE0Kj0SA7OxsajQYVFRWsAzwyMsJmWikBxquzzq/wZMKLAMLxsOnpaYyPj6OgoABFRUWifaCUSiULIVnregYGBuBwOFasrsSWuNBjWa1WdHd3hwz/XgtCmympmQcI5S9utxtms5lZtccq00gGktkppecWCssJIVhYWFjWWafVtVgjcHzCo/IlqUCShCccDwsEAqwRsH//foyOjoomFAbCa1q4XC50dHSgqKgIO3fuXJFsxXyq0+UxFVGLYVBKjyslwhNCo9GElL9QmQaVv+Tk5Kx5k262Cg9YuUvLN0GtrKxknXWLxYLJyUk2Ahfr1oKwwouDU0rUkBzh8a3XOY6DzWZDd3d3UAB1NHtuq2Gtqow2CMLJpRUThBCMj48jPT0dbW1tUX8ApU5wq2E1+cvExAQAsAqFdin5WG+yFDHg9/vD2tsVdtapCSrNpuV/PxITVJnwwoRQWzc+Po7Z2Vk0NTUFlcViE95KwuNAIIDBwUEsLCyI0hyJBA6HA1NTUygoKEBDw7IMk4ggJDwpyQQixUryF7qvSse08vLyWDW8nmUp0Z47GrIVmqDSaEYqLdJoNEHJZCv9fvIe3hoQNiZ8Ph86Ozuh1Wqxf//+ZX+8eFR4wuNR77z8/Hzs3bs3oR/eubk5DA8PY8uWLaIMXVPCk5IAVCwI5S/CTXqv1wu9Xo/CwkIm10gEklnhibV/KIxmdLlcy0YLKQHyrcWEe3hyhceDcDyMbqhu376dvdFCKJVKeDwe0a5BuKSlS9i6ujqmd0oEhBXl9PS0KEvR9bykjRTCMa3z58/D4/Ewnzoq0hVb/iLEeqzw1kJaWhrS0tKwdetWZoLKtxajHWCv18ve2zi5HUeNpBIebUzQ6mN4eBhWqxX79u1btWUeryVtIBDA0NAQbDYbs39PFNxuNzo6OpCXl8cqSrE6vvEKA5I6OI6DUqlERUUFVCoV/H5/kPyFL+KNVf4ihBSbFmKCb4JKR+BoB9jpdOKdd97BH/7wB5hMJrS2tkZ8/DNnzuDaa6/tB6AE8Cgh5PuC81cA+DWAnPde8zVCyOm1jpsUwhMuYd1uNzo7O5Gfn4/W1tY1Pyhi+tcB74+onT9/Hvn5+RFr3GKFxWJBT08Pdu7cyaQZgHgiZv7vshmJj0KpVAbtUbndbhbTaLfbodVqg9yfY0EyZSnJGC1TKBTMBFWv12Pfvn0wmUx49NFH8frrr+PJJ5/EFVdcge985ztrXpvf78fnPvc5YCl7VgfgHMdxpwgh/JnNbwD4PSHk5xzH7caSA3LVWteZcMITauvoSEwkHVAx/euAJY3b/Pw89u7dG9OIGh/hEAvtws7NzYU0PRBrKUqPE01m7nrHan8HjUbDbJr4SzQqf6FB3dFo1JJd4SVzlpZW1tdeey3+8pe/4B//8R9x4MABnDt3LqzrOnv2LGprazE8PDzy3vGexlIOLZ/wCADqSpENYDqca0sY4YXS1tEPVqQdULGWtHRMy2KxICMjQzSyC2cJ6fP50NXVBbVavWImrthL2s2IcIkn1BLNZrPBbDYzjRqd/gglfxFCCsJjKYDKUrKzs3H11VeH9TNTU1MoL+cHkkGHpaAePh4E8N8cx/0DgHQAYR08IYQn1NY5HA50dnaitLQUdXV1ET8JxSA8umdGXVbeeuutmI7HByWqlT50CwsL6OjoCNIWhoKYRGUymTAyMsJkG36/f1Msb6P9HflB3UDwkD6Vv1ACDOVSvFkrvAQllgHA7QCeIIQ8xHHcpQCe5DiugRCyaoUQd8Lz+XxwOBysgpuamsLExMSamayrIVbCM5lM6OvrW7ZnJhZWq8xooE84v78YS1A6TuTxeNDa2opAIACLxQKPx4OzZ88iPT0d+fn5LF5QRmgIh/SpRIO6FNMOpRRiGqVUXUaTZ1FaWorJyUn+l0Jlzn4KwEcBgBDyJsdxqQAKAOhXO3bcCI8uYW02G/r7+7Fnzx50d3dDpVJFFWjDR7SERzvBZrN5zU5wLAhFeNRKanFxMezfP9amBc23AIDdu3cjJSUFhBBs2bIFk5OTaG1tZZ21RMs2Eol4VFpCiYbwfaSGCGq1OuHuz8k0ABU2TKKRpbS1tWFwcBAcx1VjieiOAPi44GUTAK4C8ATHcXUAUgEY1jp2XP8SXq8XKpUKbrcb586dQ3V1NUpKVkxeCxvREB7tBGdlZUVtFBouhJXZ4uIi2tvbl1lJrYVYlrR2ux2dnZ2oqamBXh/6oSecraSyDbr8TUlJYVWLGOEyGxWh3scLFy7AbrdjenqaTYfEK6VMiGRqAIVuxwsLCxFXeCqVCj/5yU9w3XXX/QFLkpPHCSHdHMf9fwDOE0JOAfgSgF9yHPf/YKmBcTcJ42aJG+FRHZlOp8PCwgIuu+wy0YwAIyUCs9mM3t7eqJxGogG/MqPL52hEzNES3uzsLEZGRlhim8FgCOs4QtnG4uIim6t0Op2sa5mXl5fwcOn1BKVSCZVKherqaqjV6qCUMpvNBq1WywgwXuHnUtEAut3uqEYyDx48CELIDv7XCCEP8P67B8BlkR43rkvaCxcuIC0tjVnUiIVw/5g0n9VkMsV1CSsElc2MjIzAaDRGfe5Iu7SEEAwMDLBJDUpKsYT2bN26lS3baNdSp9OBEMLIL5LB8s0CftOCn1JG5S8Wi4VNKGRlZTECXO8PEiknlgFxrvDq6uqQmpqKN954I16nWRHULDMzMzPuS9hQ6O3tRXZ2dkznjoSo+OagwtlfMZ72fNPO6urqZZkVqampyMvLQ35+/rrxrIsnVmoc8OUvNPycuj/TBwnf/VlKZBEO+IQnRRVAXPfwtFqtqBMR4YJOLkSyhKXL0Fg/YHa7HQaDAVVVVdi2bVtMxwqX8Gge7Urzx2LazlPwh/ZpsJLJZFrmWbdZIxvDvdn58hcafm6xWGA0GjE0NAS1Ws3ex/UQfi7lxDJAAuYBYoIQgtHRURgMhojjGtfSzoWDqakpjI+Po7CwUBTfvHCWtNPT0xgbG0NTU1PShrT5zrpUtDs/Pw+TyYTx8XEoFAqkp6cz6y8p3QDxQrS/p0qlCpK/0H1U6lCSCIv2WMAnvGSnp4VCQict4uEGTN9QGtdIzTIjfaNjmc8NBALo7e2Fz+fD/v37MTw8HPcJCTqpQhPTVqui+JMfiZi4EHrWeTweTE1NwWw24+zZs8s0axsVYnzehfuofIt2r9cb5P6sUqmSPlHDJzyn0ykpLzwgzoRH/+DUqULMpQ0lKIVCwZawq1lKrYVo53NdLhfa29tRUlKCioqKuLicCOHxeNDe3o68vLxV7ebXOk6ioFarkZ+fD7fbjZ07dy7TrPFHtjaK9i8eCCV/4WfU0gxb6vyTjOqKf59LzfwTSFCFp1Qq4fP5RCc8n8+HyclJ6PX6iJewQkRDUitZv4tFeKGOMz8/j66uroj3J5P95KcIddOG2rPKz89PaGD3eoQwo9bj8cBoNMLr9eL8+fOskRRP+YsQfr+fVe2bnvDERmdnJzIyMqJawgoRyZKWTmxYLJaQxgfxqvDoWF6kYT4rEZ4U9tOUSmVQYpnQsVjW/oUPWkkbDAY0Njay95LKXzIzM9n+X7ziCqTsdgwkeEkrFqxWKywWC2pra1FVVSXKMcNd0tK9Qip3CUUYYlZ4hBAEAgEWq9fW1hZxpRyK8KRqCsp3LA4EArDb7Uz7B4Atf2XtX2hQwhFGNPLfy6mpKQQCATZGKOZWgpTzLIAEVnhi2TnRYJ/CwsJlIdixIBySosvJtfYKxfSx8/l8OH/+PAoLC6NyllnPoHtSK2n/0tLSRDPs3ChYae9O+F7SWV+j0Yjh4WGoVCr2XsYif5FyYhmQIMJTqVQxL2m9Xi+6urqg0Wiwf/9+DAwMiG7zvtLxCCHQ6XTQ6XRhLSfFMii12Wyw2+1oaWmJydUllOuKlPb1woVQ+yc07OQbH2xG7R8QvlPKSuHn/IAeSoCRyF82NeHRp0SsFR6trLZt24YtW7aIckwhVqrw/H5/UAh4OKW/QqFgRqfRYnJyEjqdDlqtVhQLq41AeHyEMuykeRW0Y7mZfP8oou3OrhR+TrdSwhWSy0taRL+HRwjBxMQEpqenl1VWYudahDqe0+lEe3s7ysrKUFZWFvZNE8seXiAQQE9PDwKBANra2nD27NmojsPHeie3cEAJjnYsacVCff8Sqf1L5nsthhxltfDz8fFxcBy34l7qpq7wKJRKZcQVj9frRXd3N1JSUkJWVmLnWghJijrbNjQ0IDs7O6ZjhQtqI1VcXIzKykrRqpLVmhYbFbRi0el0IX3/4qn922huxyuFn/PnqOn76fP52PkdDgdbkUkFCSO8xcXFsF9vs9nQ1dW1qn9evMK4CSEYHByE3W6POGuDf6xICY+Kp+ORhbsZCY8PofaPv2E/NDTEbO/z8vJE0f5tNMITQhh+Tt1fqJSor68P/f39LCsmEpw5cwZf/OIX4ff7MTw8/DVhPCMAcBx3K5YyLQiAdkKI0Bx0RSRMlhJO04IQgsnJSUxNTaGpqWnV9b9SqYTb7RbtWhUKBVwuF86fP4/c3NxljiORIBIyob/z9PR0zOLpta6HT8SbYU9rpb+BcMOer/1zuVzIzMxEfn5+VGllQPIjGhM9rULlL6WlpTh79izKyspw/PhxvPTSSzhz5gz++te/4iMf+Qj+9m//dtXj0HjGP/7xjygrK4NGo7ldGM/Icdx2AF8HcBkhxMJxXESjVZKRpdAUL2oBv9YfTewKb3FxEZOTk2hoaIjZJDTcCo/fEGlra4vbB5USnsPhAADWddvoFV64lVYo7Z/JZGK5Crm5ucjPzw/brXijV3irgeM4ZGVl4YEHHoDBYMCdd94Jt9uNixcvrkl4NJ6xpqaGfilUPONnAPyUEGIBAELIqhkWQkhClkLtyNdK8eJDTG0fbYwUFxeL4ogcDuHRGdytW7eivLw87jeIxWLB5OQkm3qhy7rCwkJZwMsDX68GvL9fNT09DbvdHpb2L9kBPlKR5DgcDhQWFmLPnj249tpr13x9mPGMOwCA47jXsWT//iAh5Ey415RUWQpf30btyMOFGITn8/lYsNDOnTthsVhiOh7FWoRHbd8jCR+PFoFAAHNzcyy1jH7twoULMBqNGB8fR2pqKksu20gCXjEqrWi0f8nMlEh2Yhkf0SSWhQEVgO0ArsBSmtlrHMftIYRYw/3huILjuJCyFJ/Ph56eHigUirD1bXzEKktZWFhAZ2cnKioqUFpaCovFIprMZSXCo5Mic3NzCbGc93q9aG9vh1KpRGlpKavulEolNBoNampqkJaWBqfTycw7PR7Phkkui4cl2WraP6VSidzc3KQGHiVjD49CmGcRaWJZmPGMOgBvE0K8AEY5jhvAEgGeC+ccSTEPoEvYyspKlJaWRnXMWGQpNOSGnw0rpswl1GSD3+9ne5RimB2sBRr2vW3bNvh8vmWyIH5jhW/eSZPL6CY+HUinHUwZ72Ml7d/09DTm5+fR3d3NnF/iNawvhJQaJpHq8Gg84+joKOWFUPGM/4WlEO5fcRxXgKUl7ki450gI4fErHp1Oh4mJiYiXsEJEs6SlpplOpzMo5EZ4jbFCaKlOBczl5eUoKysT5RyrQa/XY2hoiBH69PR02LIUYXKZ0HGDVn+5ubmSr/4S3Tyg2r+MjAxMTk6ivLwcJpMJXV1d8Pv9TKsWz6wKKRGex+OJSORN4xk/8pGP0Hv79yHiGf8A4BqO43oA+AF8mRBiCvscYV9NDKA3Fw2FjjWIG4ic8BYXF9HR0YGCgoKQppliTm7wydNoNKK/vx/19fWimh2EArW4N5lMaG1tZVVFLJo7YQeTWrePjo6y3Fqpetclq1tKl3ZU+1dVVRV37R//3FIhvGhw8OBBHDx4kP7vPwPL4hkJgPvf+ydiJGQPz263w+FwoLKyUrQKJxLCo7m0u3btYpWLEGIvafkxja2trTGNM4Vz49Ils1qtxr59+4I+9GIJj4WKe5q3QPVrWVlZTL8mlU5hMggv1N9rLe0f9f2LVvtHIRXCk+r8ctw/lVNTUxgdHUVaWpqoy7lwlqCEEIyNjUGv16/ZJBBzSRsIBOBwOOB2u2OOiAzHt45KXOjM70rHWOtrkYKft0DjBvnBPXRpnKwB8mTpDMMhHWHlTDN/o9X+RXLueEHqiWVAAgiPzsKePXtW1DdgreP4fD50dnYiNTU1rCaBWEtah8PBOqN1dXUxH28tYqIjaatJXBIxRsaPGwTe38AfGxuDw+FAWloa/H4/vF5vwpyLk3XDRXpe4XvH1/7ZbDZotVq2dbBWZ1/YKU0k+IRH1QBSQ9wJr7i4mL0RYgf5rATaBV5tFlcIMZa0NOOioaEB3d3dMR2LYjWyomN44Y6kCcO540mCQruhmZkZ6HQ6dHR0gOM4dgPHM2s1mYQXS5W1kvavv79/mVWTkFSSGY0odacUIIExjXTaIt6ER3NaI+0Cx3JjEEIwMjICs9kcteHASgi11KaW716vN6yRtGSbB1C7oaysLOzatQsej4ct3+x2O8tayMvLS5h8I54QU3gs1P7xk8pGR0eXORVLZUkrRS88IIGEJ/bsqxD83AcxusDhwufzoaOjA1qtdlmzQAwIiYlGNObn54dt+b4S4SULarUaW7ZswZYtW0AIYbOrnZ2dAN7fv8rKyor5QZSsLm28SEeYVCZ0Kna73dDr9QnV/lH4/X62XbFpCS9eQT4UhBAmOSkuLk5o7gMV90aydI4UfBGz3W5HR0dHRBGNwOqpZckGHTbPysoKyq2Ynp5Gf38/sxqP5gZeL3t4sUC4dfD222/D7Xajq6uLBfXk5+cjOzs77pWf1BPLgARXeGJHNSqVSiayTcRcKh9zc3MYHh4OmtaIB6iImZ6vqakp4g9Sspe0kUC4f+VwOJh4NxAIBFV/UjU9SBbRchwHpVKJqqqqIO0fNbOl2r/8/Py45NTKS1oexF7SEkLgdrsxNjYWs84t0vMODQ1hfn4+SNwbT4yNjYWcDgkX64nw+OBbjVPjTmFqGZW+hPr7b8QlbSRYSfs3PDwsqvaPgt8h3rRNi3gsab1eL5vaaGhoEJXsVrtJ6HkzMzOxb9++uN9MPp8P8/PzyM/Pj+l864HcwoFKpUJhYSEKCwtZ99JkMgXZtvOXb5thSRsJ1tL+0b3BaPdO+VKUTUt4FGItaan9+7Zt2zA3NyfqjbyayJdKXWpqasL26Y/lg0/nb9PS0mLOt1ivFd5q4HcvKyoq4Pf7YbFY2PItNTUVWq1W9G2UcJBMP7xwEUr7R40P+vv7mfYvkphGfmUrL2lFIDydTofJyUlm/24wGETPpg21HKHuKpFIXaicJBrxJfXLq6+vh06niwsxrXfCE0KpVAYt35xOJ6ampmCz2XDu3LmEWl4lyw8vlr9nSkoKiouLUVxcHFQ9U/lTdnY28vPzV33/hBVevGfHo0FCl7TRZlD4/X709vYiEAgEeefFK8iHSloIIRgYGIDD4Yh4/ywawqP5FjMzM2wUbnp6OuYJkI1GbuFAq9WiqKgIgUAAtbW1CbW8StYenlhL6VDVM9X+jYyMLNP+0XMK9/AS4QwUKSTftHA6nejo6AhphR7PMG6Px4OOjg7k5OSgpaUl4g9SKE+81UDzaAkhaG1tZR8cMcgq1FJ9M5Ag/Z1DWV6ZTKa4WV5ttGbJWto/KhwXRjRu6j28cJPL+KCjWitZK8WL8Gw2Gzo7O7F9+3YWRRftscKB2+1Ge3s7ioqKlu3XiUl4a31to2Gl348aWZSVlQW5FotleZUswkuU27FQ+2e322E2m7G4uIh33nkHL730EmZmZiJuJvIjGj/96U/ja1/7WsjXcRz3MQDHAbQRQs5Hcg5JVnh86cdqo1piE55SqcTc3BxmZ2fR3Nwc06ar0AR0JVBy3blzJ9t/iuY44VxLb28vrFYr8vPz4fF4NlR+xUpYi3iErsViWF4la0mbjPPyheMGgwHNzc0YGhrCK6+8gi984QvYsWMHbrrpJhw9enTV4wgjGtva2vD1r399Nz+i8b3zZQL4IoC3o7leyclS6FIyOzt7TSmGmIRHzS1dLpcoo2nhVHgzMzMYHR1dlVwjXRqHAtWvVVVVoaKiAvPz8xgfH4fJZILRaGTLvUS5mCQK0VRakVherXTsjbakjQQqlQpHjhzBmTNn8Ktf/Qqpqano7+9f8+eEEY1HjhxBR0eHMKIRAL4N4F8AfDmq64vmhyIF3UdZa0lrtVrR3d0d9uiUWGHcdD5VpVJh27ZtoszhrkZ4hBAMDg5iYWFhTXKNdelJw4rS0tJQXl4Ot9uNvLw8uFwuqNVqpKenw2QyMV0jvaHj6WKSKMRKPGtZXvGrP/7DYjNVeBTCzyhNLKuoqMC2bdvW/HlhRON7DY+gwBuO4/YCKCeEvMhxnHQJD1i9GqPdyenpabS0tITdORPDtHN+fh5dXV3YsWNHQpLLqNlARkZGWM2QWCo8ai+/c+dOjI6Osg40HTkqLi5GWloaKioq2CQDfzM6KysLBQUFknIwTiaEe1e0+pucnAyyvEqWLCWZXnhCso00sWwtcBynAPCvAO6O5ThJJzy/34/u7m4oFIqwrI7COWa4mJqawsTEBCPZ+fl5UQkv1FOvvb09IrOBaPfwxsfHMTs7yxxcNBoNzp8/j4yMDCwsLKCwsBDFxcUIBALs+BzHoaCggDVq7HY7y67ldzqlmF8RCvFcWnIcFxTYTS2vJiYmYDKZ4Pf7UVRUlFDLK6l44QGRd2mFEY06nQ4IjmjMBNAA4JX3/qZbAJziOO7GSBoXCVvShvrgORwOdHR0RJ3mFS3h8a2k2traWPUiZpCPkKhotbVnzx5kZWVFdJxI5S19fX3w+XxobW0FIQSEENTX1zPCzczMhNlshtVqZWLdtLQ0BAIBEELYe5qWloaqqirmYmIymdhmfnZ2NgoKCiSdXZvIvTS+5VVHRwdKSkrYdgIgnuXVapCKFx6wtJKJZE9YGNH49NNPA8Ap+n1CyDwA1tXjOO4VAP8o2S6tENT9o6GhISIC4CMawqMSkMLCwmVWUmIH+VACGR8fh16vj8rkIJIlrdfrxcWLF5Gfn4+qqir2uygUCpjNZgwMDKCxsZG5u7jdbhiNRgwPD8PpdLKbksYI+v3+oOqvqKiIjdVZrVZGgFTIS104ZABZWVkoLCxcZnnV19eHjIyMuOTVSoXwotmCEUY0Hj16FO3t7cKIxpiRcMILBAIYHByManpBiEgJjzZFVpKAiBnko1Ao4PP50NXVBY7jog7zCXdJS6u3mpoaFBYWwu/3s8pap9NhZmYGLS0tQYSr0WiChslpjODw8DA0Gg2r/jQaDSM/+n5nZmYiOzsbHMfB7XbDZDJhYGAAHo+HebDFM381HCSzW8o/r9DyamFhAWazWXTLK6kQHrDyqm41CCIaAQRHNAq+fkXkV5lgwgsEAjh//jzy8/Ojml4QIpIl6OTkJHQ63apNEaVSCY/HE9M1UQQCAQwNDaGyshIVFRVRHyccwqMzjw0NDczmm763AwMDWFxcxN69e1ddegr1aC6XC0ajEb29vfB4PMjLy0NBQQFbkvGrP5VKhZKSEmzduhWEEFgsFhgMBgwNDTELp2RYt0vRLYXjOJZXG43l1WqQEuFJVdSesD08i8UCl8uF7du3o7i4WJTjhiN1oSNbwjncUBBrSTs/P4+pqSls3bo1JrKj1+T1elf8/sTEBJu9TUlJYWTn9/vR2dmJrKws7NmzJ+Ibn8pYaI6C2WzG3Nwcc9IoLCxEbm5uUPVHCZCGzHAcx4bQp6amsLi4iJSUlCDijCfWgx4uUssrsc4rNviE5/F4JKvpTAjh6XQ6TExMICsrS1QHhbWWoIuLi2hvb8eWLVtQUVERlgQk1iXt9PQ0xsfHUV5eHratzmpYqWkRCARYitW+ffvY1ziOY5b3lZWVYVtZrQalUhl0UzocDhiNRnR3dyMQCCA/Px8FBQXIyMhgTQ/64FCr1SgtLUVWVhamp6fZv+leFh3i3wjhPRTREm04llf0/Qq1V8rPlEg01oPbMZAgwsvLy0NxcTE6OjpEHwVb6Xg0r7Wuro4t08I5XrSERwhBf38/XC4X2traMDMzI8p+YCjC83q9aG9vR25uLnbu3BnUnLBarejt7UVdXV1c7Hn4LsRVVVXMR42fQMbX7tF9P5fLBY7jWMXCcRwWFhaCwnvoRn5mZqYolZkUl7SRgG95RQhhpgf8vVK+5ZVUKjypGgcACSK89PR05pUlpiFjqA8VX8RMLZbCRbRLWuqEnJ2djebmZnAct+ZSNJJr4hOn0+nExYsXUV1djaKioqDmxOzsLCYmJtDc3JywbqnQR81ms8FoNGJiYgIKhQIFBQXgOA4zMzOor69ny21gqWlSXl7Oqhkq4qWiZ1rNRCt6TqbzsNjn5TgOWq0WWq2WbTMILa8CgUDCog6ECAQC7O+06Ss8inhHNfr9fvT0LI3eRSpiBqJb0tLksm3btgXtTYrV8eVXeGazGb29vSGbE8PDw7Db7di7d2/SpiL4Ytxt27ZhcXERAwMDsFgsSElJgU6nY9o9oewFABM9cxzHphgocYYzw7qZEMryqq+vjwWei215tRZ8Ph8rLmTCoyeLwiIqXLhcLrS3t4f0zQsXkS5pqX0VX9tGITbh6XQ66HQ67N27lz3NaQe3u7sbqampaGpqkgwZEEIwNjYGpVKJD33oQwAQJHtRq9VsuZaamgq/3x8ketZqtcuWzTTMiLrvrnUzSzVbIh5IS0uDVqtFSUkJMjIyRLe8WgtCe/dNvaSlb3C8Kjwqy4g1qjHcJS0hBKOjozCZTCvaV4nhcgIsvXcmkwkulwutra0A3m9OUGeZkpISSbnL0g5xZmYmdu7cyf7+oWQv/f39cLvdyM3NRWFhYUjZCxU90wp6fn4eJpOJ3cz8kTc+1kO2hJigo2XxsLwK59z8TFq5woP4FR6NahwaGop4vy4UwqnK/H4/urq6kJKSwuZUoz3WWvD5fBgcHATHcWhsbAxqTtjtduYsE25TJhFwu93o6OhAaWkptm7duuLrhLIX2o2kshda/aWkpCyTvWRkZATNsPLdi/nTIpsNKzUtVrO8om7GsW4XyE2LEBCzwqPEQwjB3r17RWnHhyNzuXjxIkpLS4OsbEIhVuNOmlpWWFgIl8sV1JzQ6/UYGRnBnj17JPUkdTgczCma7i2FA2E3kspeurq64Pf7meg5MzMThBDW+SWEQKlUYsuWLczBhC6bh4aGQAhh9uNiSISkjnC6tNFaXq0FmfB44JuAijHJQMmgvLwcHo9HNFX3ant4VOYS7rI5lgqPnqu+vh4qlQrvvPMOBgcHUVBQwLJEqdBYKrBYLOjv72cNlWghlL34fD4mXLbZbEGyFyq0Foqec3NzoVAoMDQ0BLfbzZK3aCUT6wiXVBGNLGUlyyvaLKLv2Vr+iELCC8fPMhlIeIUX65KW7vvQnAu9Xi/qwH+oY/EbBuHKPaIlPGpZxW9OXHLJJTCbzczhJS8vjzkVS0Gwy5fDiF1JqVSqINkLtayiVkK0MtRqtUGiZ3oDUg1oIBAIGuHSarVs708K76EYiFWHx++y19TUBFle8cN6QgnFhYRXVVUVy68SN6ybJS1tFBiNxiDXETGXycInGJ1mcLvdEctcIiU8GgnpdDrR2trKlsT035OTk8xGiy752tvbASzd9IWFhQmXbNBOrNVqTYgchp+fQG9Io9GI0dFROByOIMsqp9MJvV6Puro69vngi54dDgdMJhMb4Ocvm2N5D5M5Qyq28JhveUUfNnyhON/0QF7S8hBproUQ1HVEo9Escx2JV+eX2r7n5+dj165dcY1ppC7ImZmZaGpqYvtTCoWC7Ytt27aNLRPo8Hl1dTW76fkWT3TONZ7LNuq7BwBNTU1JWSKq1eqgzfj5+XkYjUYMDg7C7XajoqICGo0GarV6mexFrVajrKwsqGmi0+nYtAgVPUe6bZBMKUw8u9L8hw21vDKbzWxM0O12Y3Z2Fmq1OirCWyuxjOO4+wF8GoAPgAHAUULIeKS/h+SXtNQktLKyMmTXLx6EZ7fb0dnZidra2qhjGsNtWrhcLly8eBGVlZUoLi4Oak6YTCYMDg6ivr5+mc6PQnjTWywWGI1GtmzjWzyJBZ/Ph87OTuTk5KCqqkoSWjeFQoHc3Fy43W5YLBY0NDTAbrczt5jc3FwUFBQwSyuh6JlWeBzHsUpGp9Mx6/aCgoKwKujNov0TTti8/fbbWFxcxCc/+UlMT08zA9ADBw6sWfmHSiy78cYbsXv3bv7LLgBoJYQ4OY77ewA/AHBbpNedMMKjQT6RkBMV9jY0NDAZghBimnYCS2NinZ2daGxsjKksD2dJS/35du/ejaysLLaE5TgOk5OTmJ2dZXt54Z6T7ktR5w2DwYDOzk4EAgFGfrEs26gxQXl5edg29YnCxMQEjEYjW15nZ2ejrKyMVXAGg4EN4dP3glZ/fK8/Yc6HyWRioue1NGxSSA5LNOi9XV1djTNnzuATn/gEmpub8etf/xqPP/44HnvssVV/PlRi2cmTJ4MIjxDyP7wfeQvA30VzrZLU4RFCMDIyArPZvGouLT2mGBMNhBAMDw/D4/HgAx/4QMwd0LUIj7qqUFNOKhqlJgRer3dND7vVwHfeoNMKNJ9iYWEB2dnZKCwsRF5eXtjnoJblu3btikngLTbo387pdKK5uXkZ4QhlL06nk7m9+Hw+5vZCq+jVRM9CDRv92bS0NDYVsxkqvNWwuLiIm266Cf/wD/8Q1utDJZa9/faqsbOfAvBSNNcmuaYF3c9KT09fVdhLIUaFR5dodDwnnjGN/IhGfnOC/h6dnZ3Izs4OmlAQAykpKUx+wN/v4rsbFxYWrthl5VvES0n7RwPGlUplWL5//AcBP6ltenoa8/PzyMjIQEFBAfLz89meM1/2kp6ezvaxhKLnnJychPj8SQ1Ch2en07niFkys4Dju7wC0AvhwND+f0CXtWh8EOogfSapXrHt4dA+toqICpaWlePPNN0V5SociPEqsWq0Wzc3NLPOC72FXVVUlmkHqateWm5uL3NxcbN++PWTFwx/zmp6eZm7RyXLjCAX+AyLavUSVSrXMft1gMLAOOK3g0tPTg0TPAIJEz8DSFsXc3Bzm5+dZw2sz5HyESiwLN2oVCJ1YVlpauux1HMddDeAYgA8TQqIKpJZM2CgN9dmzZ09ETwelUhm1DRN1H6GaPno8MfZhhDcfndIoLy/Hli1bWHOC72G3e/fuFfcq4wmtVouKigpUVFSwPSudTof5+XkoFApwHIfm5mZJkR31BCwpKQl5c0QDvv06lb3Q5SvdBuAntQlFz1lZWUw7WV1dDbPZzLYnInEujgbJysKl5+YTHt8qKhyESiz77W9/G/QajuNaAPwCwEcJIfporzXphEcIwdDQEObn59Ha2hqxCFSpVGJxcTHi8/Kt0fnLOLq0FFNTRsO+6+rqkJ2dHdScmJmZweTkJFpaWiQx/kSFvoWFhejp6YHf70d6ejouXrzI7Mip0DdZoE7WNLAoXlCr1UFTCHQbYHR0FCqVKsjthRKfzWYDx3FsC4Gf80Gdi2PJrVgJyWyWUK9LQLzEsve8E/mJZT8EkAHgmfeIfYIQcmPE54r46qIE/+lDl3HUODMzMxP79u2L2tIpkiUt3fPx+/1obW1dtmEvZjYtAMzMzGBsbIwtB/lP4qGhISwsLCTVwy4U6D5qfn4+s8avra3F4uIiDAYDE2Pn5eWhsLAwblVLKNDGSbwcnVcCx3FsBpW+F1Tz53K5WNfWZDIxXSL/c0T96fg5H/SBwh95i7ZKSybhCc8dj8QyQsjVMV4mgCRUePSD4HQ60dnZiZqamphyFyJpWng8Hly8eBFFRUWorKwM+UcRS+ZCCMHi4iKmp6dDTk50dXVBq9VKysMOeL96CrWXmJqauizYZ2ZmhuVT8B1O4gG69N+zZ0/SlfypqakoKytDWVkZix6dnp5GSkoKent7WSUcSvSckpKC0tJS9rN8AS/N+cjPz4/ofUx2hUcf2DT4XapIOOGpVCrMzMxgYmIiZq0bPV44FZnNZkNnZ+eKmbQUYtg60c10AGhpaQlqTlD7pLKyslXtk5IBu93Olt5rVU/CYB8643rhwgVm7V5YWCia4aTBYMDIyIhklv586HQ6OJ1OXHbZZVAqlWxsraenh5kW8PWPfNkLISRI9ExzPjo6OgCAkd9aw/tSqfA8Ho+k9nqFSCjhUQ3UzMxMzCHcFOFUZLOzsxgZGUFzc/OakopYCY82J+jMq8/nY80Jm82G7u5uyenYADBLpWhkJ8IZV7fbzY7ndDrZDR3tuNv09DSmpqZEswETC1Qv6nA4gsbr+MljVPZCK+H09HRWCfNlL/QzTKto/s/S4f3Vcj6EndJEYr2YfwIJJDyv14t3330XSqUSO3fuFO2Du9oeHm2I2Gw27N+/P6x9slj28GhzYteuXcjJyYHNZsO5c+fY/o7BYEBTU1NSN/xDYWpqCtPT0xFNdawGjUaD0tJSlJaWsnE3OjVDM23pcm81CM0JknVDhwI1e/D7/avq/0LJXvjGD3zZC7Bc9FxYWBgy54OfaaHVauXEsjCRMMIzm82oqKiA0WhMSFQj3XjPyMjA3r17w15WRbuHx68iaddux44d8Pv96Ovrw9zcHFQqFYaHh9kNn+xGBZ1QcDgccSMU4bibw+EI0rnRake4ZKOE4vP5kmZOsBIIIejt7YVKpUJdXV3Yny2+7IUO4FMCs9vtyMrKCoq4FIqetVotk8y43W6YTCZm206X+cmo9NZLJi2QQMLbsmULfD4fLBZL3AmPGoRWVVVFPO8Zja3TyMgILBYL6/rS5gQhBH19fVCpVPjgBz/IBtMNBgN7StN9sESLUwOBAHp6epCSkoLGxsaENE745p78SQVq75STk8N0bn19fUhNTcXu3bsl1dShDSf6O8RybSkpKUH2S1T2EmpsTej1p1AoUFxczBp+k5OTMBgMePfdd6FWqxMqeuYHgMuE9x74FlFi5loICY8G+qxmOLDW8cIlPGozr1ar0dLSwpT4fMlNUVERKioq2M/QvS4aY2gwGNDb2wuv17tswiFeoNdWWFgYdG2JBl/nFggEYLVaodfr2ZgfdT+RSpPC7/cHyXXEBF/2AiztBZtMpiDbL+r2Qh/K/M+9RqNBXl4eampqQgZ205yPeFTK8pJ2FYht50SJgRCC8fFxzM3NBRmERopwl7RutxsXL15ESUkJysvLGYkrFAosLCygq6sLtbW1q3aE+TIPOuEwOTkJu90e1XB/OHC5XGx8L1rrq3hAoVAgIyMDQ0ND2L17NzIzM9m4m9/vZ9VOsmZVfT4fm+xIRHc9NTU1aB/UarXCYDBgaGiIzT5T2y+Px4PZ2Vls2bIFXq8XSqUySPTMz/mIh+hZblqsdsIoTUBXAyEEXV1d4DgObW1tMT3FwlnS8iUuubm5rBPLcRz7YEWa7cC3MqfD/QaDAcPDw0hNTWX7frF8SGmXOFkjbKvB6XSio6MjKACIDvjTvS76MMjKymIPg0Tsg1Iz2IqKirjPOYeCMHaRzj739vbC4/HA6/UyR5dQXn8054PjOBaPSUXPfNfiaO8bucILAX42rRhBPhRutxtOpxNlZWVsKiAWrHV9dOa3qakJaWlpQWNiExMT0Ov1MXc7+cP9ANhGf2dnJwghUVm6U/Jsbm6W3DA71f/V19cjKytr2fdD7XUZDAaWS0sfBvH4vagQe9u2batW64kEnX0uLS3FhQsXkJOTg0AggHPnzjHZCxUuC0XPKpUqSPRMcz76+/uRnp4eVc6HTHirQMwlLZWBaDQaVFZWinLMlZa0/PDtlZoTfr8fe/fuFX2fhO9rx7d0pyNNhYWFq+7PTE5OYm5uTnJJZ8D7tlPhynWEe120YqHVDt0Hpc7GscDlcqG9vZ1V8lICXWJTt2sAQRGX1PSVbgVQEhKKnoU5HzQekwqi8/Pz1zSMFRKe1Ixh+UjKklaMpgXfQPPixYuiGS+GWtLSzpxKpcLevXuDmhPU8ik3NzchdudCS3ez2Yy5uTn09/cjIyMDhYWF7OlOvffcbndciDhWzM3NBZmgRgN+oDcV6k5NTaG3t5dFOkY6pgWA7cNSN2opwefzsXxkPrkIIy6FWwH894MWHvzmh0ajYaJnOjrIz/ngx2PywSc8p9Mp7+EBwUvaWCo8frpXW1sbVCoVIykxNveFXVranNiyZQtTvwNLxEjngZPVAKAjXNTJly95USgU8Pl8yMrKQkNDg6SkHcBS1UmX/2LtwwlFvjabDUajkUmA6Hu11g1ps9nQ09MTc8ZuPOD1epl/41r7icKtAP77wf/sUNkL1fzR4oFunVA5Ff9n6dI3PT19WdNCau8ZH+uqwqNyiuzsbDQ3Ny8jUTEIj7+ktdvt6OjowI4dO5CXlxfUnLBYLOjr61tx3ynR4I93lZeX4+LFi0hPT8fi4iLOnj2bMMnLWuDbsbe0tMQ1ZYtmrFIJEHU3oaE+obYCaKA43aOVErxeLy5cuICqqqqIH7DC94MKl+nfgq+BpAUEnwDT0tJQVVUV5BNIcz7cbjfMZjMyMjLilVimAfAbAPsAmADcRggZi+gk7yGhhBdNkA+Fw+FgHmhCdxUxLZ3oH1uv17PZUjq6Q8mO7wAsFY0YBe128mMdEyV5WQs02lGhUIRlxy4m+O4mNNSHbgWkp6ezSibWJXa8QJ1+qqurRfEA1Gg0QVsjVLoyPDwMtVrNqjtqacYXPQtzPt5++21YrVbcfPPNzCh169atqK2tXfM6wkws+xQACyGkluO4IwD+BVEklgHrRJZC06ZWckMWM7mMDvm7XC7s27ePubHQm3NwcBBOpxP79u2T1GwnsNTE6enpWVZ18iUvVJcltuRlLVAHGZoHkcwqUxjqs7CwgLGxMRgMBmRkZGB6ejopweYrwePx4MKFC3HrFAtlL7QRRPNmqfkDXR0IRc9KpRI1NTV45ZVXcNtttyErKwtf+tKXsGXLFvziF79Y9dzhJJYBOATgwff++ziAn3Acx5EofKiS0qUNd0lLh8eNRuOqbshidX4DgQCGh4fh9Xpx4MCBoOYEbVykp6cnbBQrEuj1eoyOjq4pO+E4TnTJy1qgduxbtmxBWVmZKMcUC3Qo3+Px4PLLL0cgEFgWbF5QUIC8vLykNH3oHnJtbS3TJ8Yb/EYQbV7QapifdZySkoKRkRFkZGSw6m9ubg533XUXHnjggbB88cJMLCsFMAkAhBAfx3HzAPIBGCP93RK+pA33JvL7/eju7oZKpVozvUwMwqNLBrpXR5fICoVC0lmswJJdvcFgiMo+aSXJC7V1Wkvyshaojk1qkx0U4+PjMJvNaG5uZhV7qGBzOqUQj2DzlUCtxugecjIg9D3ky15cLhcLMlIoFHj22WexuLjI9vCkVhQAEsi0CAX6hy4tLQ1i/5UQK+HRtLTt27cjKysLOp0OIyMjKC4uht/vR29vb8ItxcMB7Vh7vV5RGgBCyQvNYQgleQkHDocjaCJFSljJy46PtYLNaSMolmDzlUDvASm9d3zZi0KhgMlkwpYtW/Af//EfePrpp+H1evHd7343os9hmIllUwDKAeg4jlMByMZS8yJiSI7wrFYrG38K9w8dC+Hx9wdpc+KSSy5hJgQ2m42Ne0kpVZ4aF2RkZGDHjh2i33DCm31hYQF6vT5slxe6nygFO3YhwvWy4yNUsLkw0UysRhAVPFNfRalhcnISZrOZPSj27t2Lv/zlL/jnf/5n/OUvf8GNN96Il156Kay/eziJZQBOAbgLwJsAbgbw52j27wCAW+PnRDWn9/l88Pv9eOONN3DppZcu+6DpdDrodLqIJQGjo6Os6xQuqNmAXq9HU1MTa6ZQQhsdHcX8/Dzq6+tht9uh1+thtVqRmZnJKp1k+dnR2c6tW7eKFlEYCajLi9FoDJpuoJvadAkoRWkHIYTZYm3fvl2UBwU/2NxkMoUVbL4SaJedJtxJDTqdjhnZKhQKnDlzBj/84Q/x4osvRr3sPn36NO677z6WWHbs2DE88MAD+Pa3v32IEHKK47hUAE8CaAFgBnCEEDISzbmSQnhvvfUW2tra2JMwEAiwJKw9e/ZE/ISkYshwN8RpclkgEMDu3btZ8AhtTvT09ECtVmP79u1BFR0VbxoMBphMJqjValbpJErGQJeJ/CH7ZIJONxgMBthsNqhUKrbElhrZielltxrocL/BYGDB5tTaabVzUh9HqWg7hZiensbs7CyampqgVCrxpz/9Cd/+9rdx+vTpeM0Zi/4HSgrhnTt3Dk1NTVCr1fB4POjo6EBubi5qamqi+hDqdDr4/f6w5mlpdVRQUIDKysogZTm/mxjO3qHT6YRer4fRaAQhhJFfvEZraGpXQ0NDRGHlicLo6CgMBgOys7NhsVgSJnkJB/H0slsNVANpNBphs9mYq7FwheBwONDR0SHZv+3MzAymp6dZc+e1117DN77xDbz44ovxdJBZ34Tn9/vh8/lw4cIF7Ny5k2mzamtrY+rgzczMwOVyMS3PSqDNCdri54uJ7XY7uru7o66cPB4PDAYD9Ho93G438vPzUVRUJNpkA507bWxslJzYmd882b17N6uKaUfPYDAgEAgEPRAS2cFLtJfdShCuEGiwuVarZfvIUtvvBJbiC6jQXqlU4vXXX8dXv/pVvPDCC/F+PzcG4XV0dCArKwvT09OiRDXq9XrMz89j+/btK77GaDSiv78fe/bsQXp6ehDZUREu/V6s8Pv9MJlM0Ov1sNvtyMnJYRvakTY96F6j2WxGY2Nj0nMwhAgEAuju7kZqaipqa2tXJDI6jmQwGOBwOGJOMgsXyfayWw2Li4vQ6XSYmJhAamoq2/dLZLD5WtDr9Wz6RKVS4ezZs7jvvvvw/PPPh7UKihHrn/C8Xi/eeustAEBra6sodkX0Rtq1a1fI709MTGBmZgZNTU1ISUkJmpyYmJiA0WjEnj17REnsEoLvVms2m5Geno6ioqKw5B10bzMQCKCurk4yNwEFfXjl5+dHZM/FTzKzWCxstEvsEG+qAUykaDcSUB9AWrXTvdD5+fmEBJuvBYPBgLGxMTQ3NyMlJQXvvvsuPve5z+HkyZOoqqpKxCWsb8KjIzJOpxPbt28XZSYQWNrbmpqaQn19fdDX6eymz+dDfX19UHOCetgRQhJGJnx5B39JE6qbR/ecsrOzkz6KFQq0ciorK4tJjE3fE9r1pZKXgoKCmOIsabdTqtIO6j4dygeQH2xuNBqDnE0StR1gNBpZ8HlKSgo6Ojpw77334sSJE2HNyIqE9U14dOPW5XIhMzNTtCWG3W7H6OgoGhsb2ddoAyIvLw9VVVVBzQlhZZIsMnG5XDAYDDAYDPD7/SgoKEBRURFUKpWkJzuoTiwenWLqamIwGODxeNi0RySGngsLC+js7JRst3N+fh69vb1hy3ZosLnRaBQl2HwtUBeV5uZmqNVq9PT04OjRo/j973+/4ioqTljfhBcIBOD1ejE2NsZU/WLA6XRiYGAAzc3NAN53VqHD1vz9OirrqKmpkdSok9frhdFoxPT0NKxWKwoLC1FeXo6cnBxJVXd0GZaIXAy6F0olLzTLghpYhgIVPDc2NkrSiNJqtaKvry9qjaJwOyCSYPNwYDabMTg4iJaWFqjVavT19eHuu+/GU089tWwFlQCI/sFPyg54JAYC4YDvlsKPaczIyAgiO7PZjP7+fkm2/lNSUqDRaOD1etHW1ga3242ZmRn09fUhKysLRUVFCbdzEoLasSeKTJRKZZChJ82yGBkZgUajWaaBpF52UsztAIKvL9pO+1rB5lQELgw2D/f6BgYGGNkNDw/j7rvvxpNPPpkMsosLElrhEULg8XgwMzODxcVFVFdXi3JcaoxYUlKC6enpZc0JjuOg0+kwMzODxsbGpGvCQmFmZgaTk5NoamoKuj7+jW4ymZCWlibqEz1cUFmM8PqSBTrXSiUvaWlpsNvt2Lt3r+RkO8D7D4t4eu0JO+HU1DOcByWtPOn1jY+P48iRI3j88cexb9++uFxvGFjfS1pKeOHISCKBz+fDa6+9hvz8fPYk4mdcUJfb+vp6yXnY0XCg+fl57NmzZ1XZCf+JbjQamRHjajOtYoDasTc2NkouBAhYEp6Pj48jPT09KNgo3pKXcGEymTA0NITm5uaEPSz4pp5ms5lVxAUFBcseCHRPkVaeOp0Ot956K37xi1/gwIEDCbneFbAxCG8tGUkkoLbvNpsNH/7wh4OaE3zTyWinOOIJ2kXmOA47d+6M+OakG/x6vR5er5fpuMRy76COIgsLC2hoaJDcwwJY8lOj407UrDURkpdwQeMkaQMgWaAicKPRGBRsDiCogTIzM4Obb74Z/+f//B986EMfivm8R48exQsvvICioiJ0dXUt+z4hBF/84hdx+vRpaLVaPPHEE9i7dy/99vomPGCp4zQ/P4/JyUk0NDTEdCyn04mLFy+ipqYGIyMjaGtrY0tYl8uFjo4OVFRUSLLTSTvFeXl5onSKfT4f624uLCwgNzcXRUVFUXvZEULQ29sLjuOwa9cuyT0sAAQJskORcSjJC30oxCJ5CRdSITshvF4vzGYzpqenYTabUVBQwKIyP/WpT+FHP/oRrrzySlHO9dprryEjIwN33nlnSMI7ffo0Hn74YZw+fRpvv/02vvjFL/INQDdO0yJWw06z2Yze3l7U19cjMzMT09PTePfdd1FUVASNRoORkRFJetgBYIaiFRUVy/I5ooVKpWIJVbTKoS61kTq8UOupzMxMSWoAw/GyA5YsnTIzM5GZmYmamhpWEff390cteQkXc3NzmJiYYDo2KSElJYUFPF1yySVwOp145JFH8LWvfQ2lpaXo7u5GTU2NKFnPl19+OcbGxlb8/smTJ3HnnXeC4zhccsklsFqtmJmZiVuRknDC4zguqlwLPqiN1N69e6FWqxEIBNDY2AiPx4OhoSEMDQ1Bq9XCbDZDpVJJaj6Ryjp27doVN2NHYSePzm9SOdBqDi90i6C4uFhyduxAdF52FMIgH5PJxDJsw5G8hIvZ2VlMTk6yCQWpgUqzaLd9cXERvb29ePTRR9HQ0IAXX3wRx48fx5e+9KW4X0soi/epqamNQ3hA9LIUQgj6+/uxuLiI1tZWZudEP/Q6nY5lExBCmF25y+USfZg/GphMJgwODiZUI8aP56utrWUOLzTDgj/QTyvPaGIAEwG+l11dXV1Mf8dwJC+hNvjXwszMDKamptjsqdRAXVno3LjVasWtt96Kf/qnf8J1110HAPjsZz+b5KuMH9bNkpY6XmRlZaGxsZE5EFPSowPs/LzakpISlJSUsKc5jSmMdX8rGkxPT2NqaopVpcmCVqtFVVUVy7Cgjs9OpxNerxc1NTWijfyJiXh62XEch5ycHOTk5GD79u1M8tLd3c0mYMLRtk1PT2NmZgbNzc2SJDs6bkc1qjabDbfeeivuv/9+HD58OCnXFKbFu2hIypJWoVCElWhEQY0Rq6qqsGXLFlYdKhQKuN1udHR0oKSkZMUlGP9pLtzfireol9/p3Lt3r6Q6nWq1GqWlpcjIyEB3dzeqqqowPz8PnU6H7Oxs9r4kW9qRaC87rVaLyspKVFZWsgmY0dFROByOFSUvU1NTmJubCwoDkhJoE4/ueS8sLOC2227DZz/7Wdxyyy1Ju64bb7wRP/nJT3DkyBG8/fbbyM7OjmuTUXqPIQEsFgt6enqwe/duZGVlwefzBXnYdXV1YefOnWHbSwv3t+bn51noNnUyKSgoEOUJTd2TU1JSJBntCIDZsfMdivkOL4ODg0mVdvh8Ply8eJGFCyUaKSkpbKXAl7wMDAyw98XtdrOMBymSHXWNqaurQ1ZWFpxOJ44cOYKjR4/i4x//eFzPffvtt+OVV16B0WhEWVkZvvWtb8Hr9QIA/vf//t84ePAgTp8+jdraWmi1WvzqV7+K6/UkXJbi9XoRCATwxhtv4AMf+MCqr52amsLExARr6/MnJ/R6PUZGRkTzsOM7mRiNRqjVaibqjWYJSjf/CwsLE+qwGwlmZmZYhshKv6PQ4YVWy9HkNUQKGp0pxT1F+r4MDQ3BarUiIyODvS+JkLyEC5p+Rl1jXC4Xbr/9dtx888245557kn15a2H96/D4QT4rER7txDmdTjQ0NEChUAQ1J8bHx2EymeKq/Keb+waDIeKJBrp8kGoWK7C2hm0lrOTwIrZtkdS97ICl99BisaCxsRFer5e9L9TxOl6Sl3DhdruZu3hubi7cbjc+8YlP4LrrrsNnP/tZSa44BNg4hPfmm2/iwIEDy/aHqCA3IyMDtbW1QZMTQjFsovaWaEqXXq+H3+9nnc1QchfqcyZVDSAhBIODg/B4PEF27NGA7m/p9XpmW0SbQbHcTFL3sgOW8jtsNhv27Nmz7D2MxuVFbFDvye3btyMvLw8ejwd33XUXrrjiCtx3333rgeyAjUR4/CAfCpfLhYsXL6KyshIlJSVBzQmPx4POzk4UFBSgoqIiaX8w/pNcKHehM5ONjY2SWtZQCBPZxHwP/X5/kGNvtDe51L3sAASN2631wBCaP8QieQkXlOxodez1enH06FHs378fX/nKV9YL2QEbgfBorsW7776Luro6tkTkB3BnZWWF9LDbtm2bpCQT/OwKs9kMANi1axcKCgqS3tkUgnY6c3Nz427PLbzJU1NTWTNotf1QqXvZ0Y67y+WKujrmRzhGInkJF9Q5qKamBgUFBfD5fLj33ntRV1eHb37zm+uJ7ICNRHgdHR2oqalBRkYGpqeng6yH+GRHxbq0nS41EEIwNDQEp9OJkpISmEwmWK1WyXjYAeLZsUcLh8PBmkEcx6GwsBBFRUVB+6HUPkmK4d3A+39nuhUgFjlR8ltN8hLJ8S5cuIDq6moUFhbC7/fjc5/7HCoqKvDtb397vZEdsJEIr7u7G1u3boXRaITdbmd7IXyym5ycZE4YUhq+pqCCZ41GE7RE5MtdTCaT6HKXSEDt2Gtra+MVlhwRQjm8qFQqzMzMxNUrLhbQfU+fzxfzhMdKiNXlhcafVlZWMr3pfffdh9zcXPzLv/yL5FYcYWL9Ex61ee/t7YXNZmPqdmFzgr+xnuwKKRRoZkZxcfGqcXViy10iQSLt2KOBz+fD0NAQZmdnoVar2TC/VHzsgPcVA4QQ7Ny5MyFVUijfQ0p+oZb6VKtYXl6O4uJiBAIBfPnLX4ZarcaPf/xjybyXUWBjEJ7dbsdbb72FoqIi1NXVBTUnfD4fOjs7JZvWBbxfNUWzpxiL3CUSUDtxsXSK8YBOp8Pc3BxzPOFXOFTXFq7DSzxAk+0UCgV27NiRVHnJSpKXQCDAhNlUHH3s2DF4PB789Kc/FY3szpw5gy9+8Yvw+/349Kc/ja997WtB35+YmMBdd90Fq9UKv9+P73//+zh48GCsp13/hGez2XD+/Hnk5OQgOzsbxcXFyzzsKisrRbNNEht0Y12MqimU3IVq2mKBXq/H6OgompqaJGl3DryvYduzZ8+yCp7v8GIymdZ0eIkHqARKpVKJ3tGOBfxuuNVqhdfrRXFxMSorK5GamooHH3wQJpMJv/zlL0VbGfn9fuzYsQN//OMfUVZWhra2Njz11FPYvXs3e80999yDlpYW/P3f/z16enpw8ODBVW2hwsT698PTarVoaWmByWTCzMwM1Go1cnNzmc20VJdfwPtEIlZITGpqKsrLy1FeXs7kLtSOnnbvInV30el0mJ2dxd69eyVpTUQ7nU6nE42NjSErkFAOLwaDIaTDS7yusaenBxqNBtu2bZMM2QFgub15eXm4ePEiW2FcffXV8Pl80Gq1eOaZZ0TdBjp79ixqa2tRU1MDADhy5AhOnjwZRHgcx8FmswFYKgqSMQYYDhJOeIFAgO1jaTQa6PV6dHd3IxAIYOfOnZLsxAJLJbvBYIgbkaSkpLB5Ub/fD6PRuMzdJTc3d8Wbj2Zj2O12tLS0SHLfk+9l19DQEDaR8If5+Q4vi4uLok80UK1iWlqaJGMBgKVr5HsWEkJw00034d1338Xll1+Oo0ePoqioCL/97W9FOV8ozzqeKzEA4MEHH8Q111yDhx9+GA6HAy+//LIo5xYbCSe83//+93jooYfw0Y9+FDfccAN+85vf4EMf+hA+9KEPwWQyYWxsjIV0S0HSQRsobrcbLS0tCdkAViqVKC4uZhvQFosFs7OzQe4u+fn57FroXhMhRLImBbRqUqvVMe2HUYeX0tJSpoPU6XTo7e2N2eGFb0FFqxmpIRAIMAE+JbuHH34Y7e3teOaZZ5CSkoL777+fDegnCk899RTuvvtufOlLX8Kbb76JO+64A11dXZJrmCSc8P7u7/4ON9xwA06cOIEjR46wzlNlZSXzsrPZbJibm2MOJsXFxSgoKEg4+VGr8/T09IgqEjGxlrtLQUEB5ubmJBtUBLxPJJmZmaiqqhLtGoW2X7E4vFAiyc7OjrswO1rQa8zJyUF5eTkIIfjFL36BN954A8ePHw/6PcVchYTjWffYY4/hzJkzAIBLL72UyY+kNkue8KYFsCREveaaa/CZz3wGN998M06fPo0TJ06gt7cXV155JQ4fPswCeex2O5N0pKWlsa5mvDt3VKy7ms9eMkEIgdVqZcEofLcOKWkW/X4/2tvb2UhgIhAqvGc1hxe6RMzNzRUlxyEeIISwBL6qqioQQvD444/j9OnTeO655+LanPL5fNixYwf+9Kc/obS0FG1tbfjtb38bFM597bXX4rbbbsPdd9+N3t5eXHXVVZiamor14bb+u7TA+6p1YS6ty+XCH/7wB5w4cQIXLlzA5ZdfjsOHD+PSSy+FQqGAw+HA3Nwc07MVFxejsLBQ9D01OrwuFbFuKLjdbrS3t6OyshLFxcUh5S5FRUVJ7dJSrWKyvOwoQjm80HEuSnaJMheNBoQQdHd3Q6vVsqX2b37zG5w4cQKnTp1KyGTK6dOncd9998Hv9+Po0aM4duwYHnjgAbS2tuLGG29ET08PPvOZz2BhYQEcx+EHP/gBrrnmmlhPuzEILxy43W68/PLLOH78OM6dO4dLL70UN910Ey677DKkpKSwcSWDwQCVSsVu8FirG6vVit7eXjQ0NEi2gUJni3fs2BHS+DRecpdIIFUvO6HDi8/nQ3FxsaSkJ3zQvc/U1FRs27YNAPD000/jySefxAsvvCBZjaVI2DyEx4fX68Urr7yC48eP4/XXX0dbWxsOHTqEK664Amq1Gi6XC3NzczAYDFAoFIz8ItVszc3NYWxsTNL6NaoDDJeQqdxFr9czuUtRUZFoYd2hsB687Px+Py5cuMCqvFgcXuIFqgVMSUlBbW0tOI7DiRMn8Oijj+KFF16Q7ANZRGxOwuPD5/Phr3/9K5555hm89tpraGpqwqFDh3DVVVchNTUVi4uL0Ov10Ov1IIQw8lut7CeEYGJiAiaTCXv27JGkfg14P/Us2gF7KnfR6/VBYd2ryV0ixXrwsqOBUCUlJWypHa3DS7xAE/oUCgWrPk+dOoWHH34YL774omTfW5EhEx4f1Ej0xIkT+NOf/oRdu3bh8OHDuOaaa6DVauF2uxn5+f1+Rn58rzr6wfL7/airq5NcG50iHDv2SBAIBGA2m6HX61l1I5S7RArqZSfl7QA6d1paWrqqc0woh5dE2bdTvSIAJuF56aWX8MMf/hCnT58OO79lA0AmvJUQCARw/vx5HD9+HP/93/+NmpoaHDp0CB/96EeRmZnJBKt8l46CggKMjIwgKytLsnO7wJLo2Wg0orGxMS7daVrdzM3NwWw2R+XuQidlpDy76/V6cfHiRVRUVKC4uDjsn6OzrPSzQ01f47EtQBt6Pp8Pu3btAsdxePnll/Gd73wHp0+flmwTLU6QCS8cBAIBtLe34/jx43jppZdQWlqKQ4cO4eDBg8jJyYHX68XU1BRGR0ehUqmwdetWFBUViWbCKBboh39xcRH19fUJqT6jcXeRupcd8L5XXKxNFJ/Px0xf6baAmA4v1HOP2lC9+uqr+OY3v4kXX3wxIpLeIJAJL1LQlv7x48fx4osvIj8/Hx/4wAfw1FNP4dlnn0VpaWlQ145v2Z5M8gsEAmx4PZlOHdSmiDaEaMeXNnUMBgNGRkbQ3NwsSS874P2OMTXGFAtCD7tYHV6Gh4exuLjIDEZff/11fPWrX8ULL7wg2dnUOEMmvFhACMH//b//F1/5ylewfft2pKSk4MYbb8QNN9zA1Pomkwlzc3NYWFhAXl4eiouLE548lUg79kgglLukpqbC6XRi3759km300HyHbdu2xXU5SAhhIvloHF5GR0dZTgbHcXj77bdx//3349SpU6v6LW5wyIQXC0ZHR/Hxj38cJ06cQElJCUZGRnDixAmcPHkSKSkpuOGGG3Do0CGUlJSAEMKWLjabLS4dzVCgEx6lpaWSfqqPj49jamoKaWlpcLvdCZG7RAq3242LFy8mRR5DHV4MBgMIIUFxlkKMjY3BZrOxUKB33nkHn//853Hy5EnRHnhr+dkBS3PuDz74IDiOQ1NTk2jmAzFAJrxY4fF4lu1FEUIwOTmJEydO4L/+67/g9/tx/fXX4/Dhw2xm0WKxYG5uDvPz88zHT2xnXqnZsa+EsbExWK1W5mUnlLvwnYuTncm6kjg7kfB4POz9ETq8TE5OMl9AhUKB9vZ23HvvvXj22WdRW1sryvnD8bMbHBzErbfeij//+c/Izc2FXq+XgmBcJrx4gxCCmZkZPPvss3juuefgcDhw/fXX49ChQ2ysx2KxQK/Xw2KxiCLnAN6XdEg1zxZYem+Gh4fhcrlWbKLEQ+4SKRYXF3Hx4kUWQC0l8DNrTSYTOI5DSUkJSktLMTIygk996lN45plnsHPnTtHO+eabb+LBBx/EH/7wBwDA9773PQDA17/+dfaar3zlK9ixYwc+/elPi3ZeEbD+DUClDo7jsHXrVnz+85/H5z//eej1evzXf/0X/vEf/xFmsxkHDx7EoUOH2AeSyjmGhoaQkZGB4uLiiJX6VqsVfX19ko0nBN7XKxJCVnWOUSgUTPJDDQ747i7xDjOiVbJUhc/UyMDj8cDtdqOyshKnTp3CQw89BLfbjS984QuiV1bh+NlR3d9ll10Gv9+PBx98EB/96EdFvQ4pQCa8NVBUVIR77rkH99xzD8xmM06ePIkHHngAMzMzuOaaa3DTTTexrprNZoNer8fIyAi0Wm1YNzffRVmq42y0Y6xWq9mIUzjgOA65ubnIzc1lcpe5uTmMj4/HJcyIkl1dXZ1kXbMBYHp6Gnq9Hk1NTVAqlfjgBz+IRx55BD/4wQ/Q19eHgwcP4p/+6Z9www03JOyafD4fBgcH8corr0Cn0+Hyyy9nVlQbCTLhRYC8vDx88pOfxCc/+UnMz8/j+eefx/e//32Mjo7i6quvxuHDh9HU1ITa2lqmZRsfH4dGo2Gefvxu5tTUFGZmZiRrxw6878FGxdnRguM4ZGZmIjMzE7W1tUzu0t7eHlLuEinoSBsNcpcqZmZmMDMzg+bmZiiVSoyNjeGOO+7AE088gX379gEAvvGNb2CNraaIEI6fXVlZGQ4cOICUlBRUV1djx44dGBwcRFtbm2jXIQWItoe3VhfI7XbjzjvvxDvvvIP8/Hz87ne/k5TkIhbY7Xbm6dff3888/VpbW5fZWlFnl8XFRSwsLIQMsZEKqJddYWFhXKURsbq7OBwOdHR0SHqkDQBmZ2eh0+mYBf/k5CRuu+02/OIXv8CBAwfidt5w/OzOnDmDp556Cr/+9a9hNBrR0tKCixcvJtv8QZpNi3C6QD/72c/Q0dGBRx55BE8//TSee+45/O53v4vx8qUH6ul3/PhxtLe3M0+/Sy65BEqlEna7HZ2dnfB6vUyoGo2zS7yRLC87YUdzLbnLepjfBZaceCYmJtDS0gKVSoXp6WnccsstePjhh/HBD34w7udfy8+OEIIvfelLOHPmDJRKJY4dO4YjR47E/brWgDQJL5wu0Ec+8hE8+OCDuPTSS+Hz+bBlyxZmVrlR4Xa78cc//hHHjx/H+fPnceDAAfT39+PgwYP44he/GOTsIhXTTkA6Xnah5C5FRUXIyclhbthdXV3Ys2cPMjIyknada8FgMGBsbAzNzc1ISUnB7Owsbr75Zjz00EO48sork315UoY0u7ThdIH4r1GpVMjOzobJZJK03ixWaDQaXH/99bj++uthMplwzTXXIDs7G7/97W8xODiIw4cP48Mf/jAqKyuZswtNcCssLERxcXHCZ1Ol5GUnDDMym82YmZlBX18ftFotbDYbmpubJU12RqMRo6OjaGlpQUpKCgwGA2655RZ8//vfl8kuCZCbFgnC/fffj69//eu4+eab4fP58Je//AXPPPMMjh07hubmZhw+fBh/8zd/g/Lycng8Huj1evT29sLn8yXMsVjKXnZ8ucv8/Dw6OzuRl5fHQpaoHCjeWSeRwGQyYWRkhJGdyWTCLbfcgm9/+9ti2J/LiAKifDrC6QLR15SVlcHn82F+fj7pFUQi8dhjj7GbUaVS4corr8SVV14Jv9+PN954AydOnMC3vvUt7N69G4cPH8bf/u3foqysjDkWDwwMwOPxoKCgAMXFxUhPTxd1O2C97IVRG6p9+/YhLS0tSO4yNjYWF7lLNDCbzRgaGmJkZ7Vaccstt+DYsWM4ePBg0q5rs0OUPbxwukA//elP0dnZyZoWzz77LH7/+9/HePkbC4FAAOfOnWOefrW1tczQNDMzEz6fD0ajEXNzc3C5XKLNr64HLzvgfYH2ajZUa7m7JAIWiwUDAwPMQcZms+FjH/sY7rvvPtxyyy0Ju44NAGk2LYC1u0CLi4u44447cOHCBeTl5eHpp5+WbNixFEA9/Z555hm89NJLKC8vZ55+2dnZyzb08/PzUVxcHLGt1XrwsgOWSKS/vz8igXYywowoKbe0tECj0WBhYQG33HIL7r33Xnz84x+P23k3KKRLeDLiB0IIurq6mKdfYWEhDh06hOuuuw75+fnw+/0wm82Ym5uD3W5f1s1cCQaDAaOjo2hqapKcLIYPSsqxTKNEKneJBrRSptfpdDpx66234q677sJdd90lyjk2GWTC2+ygM63Hjx9nyVXU06+wsBCEkKDh/ZycHGZrxR/en52dxeTkJJNKSBUmkwlDQ0OiGowKXYvDfUCsBpvNhp6eHlYpu1wuHDlyBLfeeis+85nPiHLdmxAy4Qmx1oTHv/7rv+LRRx+FSqVCYWEhHn/8ccmmy0cKQghGRkZw/PhxnDx5EhqNhnn6bdmyhQ3vz83NwWq1IisrC8XFxXC5XGyWU0pdTSGMRiOGh4fR0tIStwaE0N0lOzsbRUVFyMvLC9vdxW63o7u7m5Gd2+3GJz7xCVx33XX47Gc/u6G1pnGGTHh8hDPh8T//8z84cOAAtFotfv7zn+OVV17ZkBMeNGqS2loFAgHccMMNOHz4MMrKygAsLbkGBgawsLDAdH5SyWAVgi63m5ubE9Zt5bu7mM3msCzbaXe7qakJWq0WHo8Hd955J6688krcd999opFdOAaeAHDixAncfPPNOHfuHFpbW0U5dxIhEx4f4Ux48HHhwgV8/vOfx+uvv56wa0wG+J5+zz77LFwuFw4ePIiBgQGUlpbi2LFjzNzAZDJBq9UycwMpkB81XUjmcjuUZTt9jygBU7Kjtl5erxdHjx7F/v378ZWvfEU0sgvnwQ4sVZrXXXcdPB4PfvKTn8iEFwLSXc+EgXAmPPh47LHHcO211ybi0pIKvqff5z73OczNzeHjH/84ZmZm0NvbC41Gg8OHD2PHjh3M2WVubg6jo6NIS0tjOrZkLHfpzGmy9xY5jkNWVhaysrKYu4ter2fuLtnZ2Zibm0NTUxPS09Ph8/lw7733orm5WVSyA4CzZ8+itraWqRqOHDmCkydPLiO8b37zm/jqV7+KH/7wh6Kde6NhXRNeJPjP//xPnD9/Hq+++mqyLyWh4DgO//qv/4q2tjZ8//vfZ55+3/jGNzA7O4uPfOQjuOmmm1BXVxdEfu+88w6ragoLCxNCPjMzM5iammID9lJCeno6qqurUV1dDYvFgs7OTmg0Gvz85z+HxWKBTqfDrl278I1vfEP0PbtwHuzvvvsuJicncd1118mEtwqk9amKEOFMeADAyy+/jH/+53/Gq6++Kmn5Rbzwla98hc0s5+fn4+jRozh69CisViuef/55fPe738XY2Biuvvpq3HTTTWhsbMS2bdtYVXPhwgVma1VUVBSXPbXp6WnmEyc1suPD5XKhv78fLS0tyMzMRGVlJe6//3709PRgdHQUAPDpT386odZngUAA999/P5544omEnXO9QrqfrDDQ1taGwcFBjI6OorS0FE8//fSypKULFy7g3nvvxZkzZ6QQSpIUrGTQkJOTgzvuuAN33HEH7HY7XnzxRfz4xz9Gf38/rrrqKhw6dAitra2orq5mnV26pKPLXjEmGKampjA7O8tMMaUKaqywe/duZGZmIhAI4Ac/+AEqKyvx7LPPwul04syZM3A6naKed60HO3WNueKKKwAsSY5uvPFGnDp1aiPs44mKdd20ANae8Lj66qvR2dmJkpISAEBFRQVOnTqV5KuWNlwuF1566SWcOHECHR0d+PCHP4xDhw4xTz9qa2UwGBAIBFjlF82khk6nC7I7lypoMBC1jw8EAjh27Bg8Hg9++tOfxjWgKJzRTT6uuOIK/OhHP9oIZCd3aWUkFtTT75lnnsE777yDD3zgA7jppptw2WWXQaVSMWcXvV7PnF2Ki4uh1WrXPPbExARMJhMaGxslTXY08pG6yAQCAXzrW9+C2WzGf/zHfyTk2td6sPMhE94qB5QJT0a48Hg8+J//+R8cP34cb775Jvbv34/Dhw/j8ssvh1qtZs4uc3Nz8Hg8bHY1lF/d+Pg4LBYLGhsbExbfGA08Hg8uXLiA7du3Iy8vD4QQfPe738XExASeeOIJSRP1BoBMeFLCJhWDAlhaZr322mt45pln8Je//AUtLS04fPgwrrzySqSmpsLr9bLZVersUlxcjIyMDIyPj2N+fp6FT0sVlOyoGSohBD/60Y/Q19eHJ598UtLNlQ0CmfCkgk0sBl0Gv9+P119/HSdOnMCf//xn1NfX4/Dhw7j66quh1WrZ7Orc3BwsFgtUKhXq6+uRnZ0t2bErr9eLCxcuoKamhmXsPvzwwzh37hyefvppSc8fbyCI/uGQ7uNV4uCLQdVqNRODCkHFoMnOqYgnlEolLr/8cvz7v/872tvbcd999+H8+fO46qqrcMcdd+DkyZNIS0vDc889B0IIamtrodPp8NZbb6G/vx9Wq1XUWMJYQcmuurqakd0vfvELvPHGG3jqqadkslvHkAkvSoQSg05NTQW9hi8G3SxQKBS45JJL8KMf/QgXLlzAsWPH0N3djf379+Pll1/G1NQUUlNT0dDQgAMHDiA/Px9TU1N466230NfXB7PZjEAgkLTr9/l8LMCIus889thjePnll/H73/8+qS7KMmKHTHhxAhWDPvTQQ8m+lKRBoVCgpaUFgUAAH/3oR/H4449jcnISN9xwAz72sY/hP//zP6FQKFBfX48DBw6gsLAQc3NzePvtt9HT0wOTyZRQ8qNkV1FRwTSbTz75JJ5//nmcOHFiQ1fpmwXyrmuUkMWg4cHv96O2thaf/vSnoVAo0NTUhG9961vo6+vD8ePH8bGPfQzZ2dm48cYbcf3112PXrl0AwGytBgYGkJWVxVxL4tXkoKHjZWVlKC4uBgA89dRT+N3vfocXXnhB0m7QMsKH3LSIEptYDCoqCCEYHh7G8ePHcerUqWWefsCSrRV1LaGWTWI6u/j9fly8eBFbt25lAvXjx4/j0UcfxYsvvijpUKMNDrlLKyVsUjFo3EAIwfj4OPP0A4Drr78+yNOPWjYZjUakpaUxy6ZoJSK0sisuLmYV+qlTp/Dwww/jxRdflFxc5SaDTHgyNgeop9+JEyfw7LPPYnFxEddffz0OHTqE6upqAEsJZXNzczAajdBoNGy+N9wuKg1KKiwsZIT60ksv4Yc//CFOnz6NvLy8uP1+MsKCTHgyNh8IIdDr9Xjuuefw7LPPwmq14uDBgzh06BB27NgBjuOYs4vBYIBKpWK2Vit1VQOBADo6OpCfn8+67S+//DK+853v4PTp0ysaLkSDzRxDECPEF2kSQlb7R4ZE8dJLL5EdO3aQbdu2ke9973shX/O73/2O1NXVkd27d5Pbb789wVcYPxiNRvLoo4+Sa6+9luzbt48cO3aMnDt3jiwsLBCHw0GMRiPp7u4mr7zyCnnttddIX18fMZvNxOFwEIfDQex2O3njjTdIT08P+9rp06fJvn37yOzsrKjX6vP5SE1NDRkeHiZut5s0NjaS7u7uoNf8+c9/Jg6HgxBCyM9+9jNy6623inoN6xhr8VPE/8iEtw4Rzk00MDBAmpubidlsJoQQMjc3l4xLjTssFgv59a9/TQ4dOkSam5vJl7/8ZfL6668Tu91OHA4HMZlMpLe3l7z66qvk1VdfJd3d3eSvf/0r6e7uZmT3hz/8gbS0tJCpqSnRr++NN94g11xzDfv/7373u+S73/3uiq9/9913yQc+8AHRr2OdQnTCk2Up6xDhWH7/8pe/xOc+9znk5uYCwIb1AszJycGdd96JO++8EzabjXn6DQwM4G/+5m9w+PBh7Nu3DxUVFXA4HLh48SIIITAajXj88cfR1NSEhx56CM8//zy2bt0q+vXJMQTSgiw8XocIZ8pjYGAAAwMDuOyyy3DJJZfgzJkzib7MhCMrKwu33347nnnmGbz55pv44Ac/iEceeQSXXnopvvzlL+PjH/843nrrLXzwgx/Evn37QAjB17/+dXAch9/85jfo7+9P6vXTGIIvf/nLSb2OjQy5wtug8Pl8GBwcxCuvvAKdTofLL78cnZ2dm0ZmodVq8bGPfQwf+9jH4HQ68b/+1/+C2WzG008/DZ1Oh8bGRvz3f/83XnnlFRQUFOD555/HCy+8gJ07d4p6HXIMgbQgE946RDg3UVlZGQ4cOICUlBRUV1djx44dGBwcRFtbW6IvN+n48Y9/jD179uAHP/gBvF4v/vznP+OBBx7AL3/5S9TW1gIA7rjjjricW44hkBjW2OSTIUF4vV5SXV1NRkZGWNOiq6sr6DUvvfQSufPOOwkhhBgMBlJWVkaMRmMyLjfpsFgsJBAIJO38L774Itm+fTupqakh3/nOdwghhHzzm98kJ0+eJIQQctVVV5GioiLS1NREmpqayA033JC0a5UYRG9ayDq8dYq1pjwIIfjSl76EM2fOQKlU4tixYzhy5EiyL1uGjEggC49lyJCxaSAbgMqQIUNGtJAJT4YMGZsGMuHJiBhnzpzBzp07UVtbi+9///vLvj8xMYErr7wSLS0taGxsxOnTp5NwlTJkLIe8hycjIoQTXnTPPfegpaUFf//3f4+enh4cPHgQY2NjybtoGesV8h6ejOQinPAijuNgs9kALJl3xmNkS4aMaCALj2VEhHBmQx988EFcc801ePjhh+FwOPDyyy8n+jJlyAgJucKTITqeeuop3H333dDpdDh9+jTuuOOOpCaRyZBBIROejIgQzljbY489hltvvRUAcOmll2JxcRFGozGh1ylDRijIhCcjIvBnQz0eD55++ull+R0VFRX405/+BADo7e3F4uIiCgsLk3G5MmQEQSY8GRFBpVLhJz/5CT7ykY+grq4Ot956K+rr6/HAAw/g1KlTAICHHnoIv/zlL9HU1ITbb78dTzzxBDhOfLfueGMt+Y3b7cZtt92G2tpaHDhwQO5ErwPIshQZMkIgHPnNz372M3R0dOCRRx7B008/jeeeew6/+93vknjVGw6yLEWGjEQgHPnNyZMncddddwEAbr75ZvzpT3/CGgWEjCRDJjwZMkIgHFdp/mtUKhWys7NhMpkSep0yIoNMeDLWHY4ePYqioiI0NDSE/D4hBF/4whdQW1uLxsZGvPvuuwm+QhlSxVp7eDJkSA4cx10OYAHAbwghy1iP47iDAP4BwEEABwD8OyHkQITnuBTAg4SQj7z3/18HAELI93iv+cN7r3mT4zgVgFkAhUS+qSQLucKTse5ACHkNgHmVlxzCEhkSQshbAHI4jiuJ8DTnAGznOK6a4zg1gCMATglecwrAXe/9980A/iyTnbQhE56MjYhSAJO8/9e997WwQQjxAfg8gD8A6AXwe0JIN8dx/x/HcVR4+BiAfI7jhgDcD+BrMV+5jLhCnqWVIWMFEEJOAzgt+NoDvP9eBHBLoq9LRvSQKzwZGxFTAMp5/1/23tdkbHLIhCdjI+IUgDu5JVwCYJ4QMpPsi5KRfMhLWhnrDhzHPQXgCgAFHMfpAPy/AFIAgBDyCJaWoQcBDAFwAvhkcq5UhtQgy1JkyJCxaSAvaWXIkLFpIBOeDBkyNg1kwpMhQ8amgUx4MmTI2DSQCU+GDBmbBjLhyZAhY9NAJjwZMmRsGvz/9TyXFF6u1DcAAAAASUVORK5CYII=\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "from mpl_toolkits.mplot3d import Axes3D\n",
    "from pylab import mpl\n",
    "import math\n",
    "import matplotlib.pyplot as plt\n",
    "\"\"\"多元线性拟合\n",
    "其中拟合数据为x1=2, 4, 5, 8, 9\n",
    "             x2=3, 5, 7, 9, 12\n",
    "             y=48, 50, 51, 55, 56\n",
    "\"\"\"\n",
    "x1 = [5.41,2.99,3.76,3.95,3.68]\n",
    "x2 = [7.2,7.3,6.3,7.4,8.3,7.8]\n",
    "y = [37.1446,28.7479,23.8610,27.9866,38.1060]\n",
    " \n",
    "\"\"\"完成拟合曲线参数计算前相应变量的计算\"\"\"\n",
    "def multivariatel_liner_fitting(data_x1,data_x2,data_y):\n",
    "    size = len(data_x1)\n",
    "    sum_x1 = 0\n",
    "    sum_x2 = 0\n",
    "    sum_square_x1 = 0\n",
    "    sum_square_x2 = 0\n",
    "    sum_x1x2 = 0\n",
    "    sum_y = 0\n",
    "    sum_x1y = 0\n",
    "    sum_x2y = 0\n",
    "    i = 0\n",
    "    while i < size:\n",
    "        sum_x1 += data_x1[i]\n",
    "        sum_x2 += data_x2[i]\n",
    "        sum_y += data_y[i]\n",
    "        sum_x1y += data_x1[i]*data_y[i]\n",
    "        sum_x2y += data_x2[i]*data_y[i]\n",
    "        sum_x1x2 += data_x1[i]*data_x2[i]\n",
    "        sum_square_x1 += data_x1[i]*data_x1[i]\n",
    "        sum_square_x2 += data_x2[i]*data_x2[i]\n",
    "        i += 1\n",
    "    return [[size, sum_x1, sum_x2, sum_y]\n",
    "        ,[sum_x1, sum_square_x1, sum_x1x2, sum_x1y]\n",
    "        ,[sum_x2, sum_x1x2, sum_square_x2,sum_x2y]]\n",
    " \n",
    " \n",
    "\"\"\"利用高斯消元法，求解拟合曲线的参数\"\"\"\n",
    "def calculate_parameter(data):\n",
    "    #i用来控制列元素，line是一行元素,j用来控制循环次数,datas用来存储局部变量。保存修改后的值\n",
    "    i = 0;\n",
    "    j = 0;\n",
    "    line_size = len(data)\n",
    " \n",
    "   #将行列式变换为上三角行列式\n",
    "    while j < line_size-1:\n",
    "        line = data[j]\n",
    "        temp = line[j]\n",
    "        templete=[]\n",
    "        for x in line:\n",
    "            x=x/temp\n",
    "            templete.append(x)\n",
    "        data[j]=templete\n",
    "        #flag标志应该进行消元的行数\n",
    "        flag = j+1\n",
    "        while flag < line_size:\n",
    "            templete1 = []\n",
    "            temp1=data[flag][j]\n",
    "            i = 0\n",
    "            for x1 in data[flag]:\n",
    "                if x1!=0:\n",
    "                   x1 = x1-(temp1*templete[i])\n",
    "                   templete1.append(x1)\n",
    "                else:\n",
    "                   templete1.append(0)\n",
    "                i += 1\n",
    "            data[flag] = templete1\n",
    "            flag +=1\n",
    "        j += 1\n",
    " \n",
    " \n",
    "    #求相应的参数值\n",
    " \n",
    "    parameters=[]\n",
    "    i=line_size-1\n",
    "    #j标识减去元素个数\n",
    "    #flag_rol标识除那一列\n",
    "    flag_j=0\n",
    "    rol_size=len(data[0])\n",
    "    flag_rol=rol_size-2\n",
    "    #获得解的个数\n",
    "    while i>=0:\n",
    "        operate_line = data[i]\n",
    "        if i==line_size-1:\n",
    "            parameter=operate_line[rol_size-1]/operate_line[flag_rol]\n",
    "            parameters.append(parameter)\n",
    "        else:\n",
    "            flag_j=(rol_size-flag_rol-2)\n",
    "            temp2=operate_line[rol_size-1]\n",
    "            #result_flag为访问已求出解的标志\n",
    "            result_flag=0\n",
    "            while flag_j>0:\n",
    "                temp2-=operate_line[flag_rol+flag_j]*parameters[result_flag]\n",
    "                result_flag+=1\n",
    "                flag_j-=1\n",
    "            parameter=temp2/operate_line[flag_rol]\n",
    "            parameters.append(parameter)\n",
    "        flag_rol-=1\n",
    "        i-=1\n",
    "    return parameters\n",
    " \n",
    "def calculate(data_x1,data_x2,parameters):\n",
    "    datay=[]\n",
    "    i = 0\n",
    "    while i < len(data_x1):\n",
    "        result = parameters[2]+parameters[1]*data_x1[i]+parameters[0]*data_x2[i]\n",
    "        i += 1\n",
    "        datay.append(result)\n",
    "    return datay\n",
    " \n",
    " \n",
    "\"\"\"完成函数的绘制\"\"\"\n",
    "def draw(data_x,old_y,new_y):\n",
    "    #创建绘图函数对象\n",
    "    fig = plt.figure()\n",
    "    #创建Axes3D对象，让其包含图像3D坐标\n",
    "    ax = Axes3D(fig)\n",
    "    ax.scatter(data_x[0], data_x[1], old_y, color='red')\n",
    "    ax.plot(data_x[0], data_x[1], new_y, color='black')\n",
    "    mpl.rcParams['font.sans-serif'] = ['SimHei']\n",
    "    mpl.rcParams['axes.unicode_minus'] = False\n",
    "    plt.title(\"多元线性拟合数据\")\n",
    "    plt.show()\n",
    " \n",
    " \n",
    "data = multivariatel_liner_fitting(x1, x2, y)\n",
    "data2 = calculate_parameter(data)\n",
    "newY = calculate(x1, x2, data2)\n",
    "draw([x1, x2], y, newY)\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "74927c2e",
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "id": "f18f5f8e",
   "metadata": {},
   "outputs": [
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "UsageError: Cell magic `%%基于BP神经网络的粮食产量预测` not found.\n"
     ]
    }
   ],
   "source": [
    "\n",
    "clc\n",
    "clear all\n",
    "close all\n",
    "format short\n",
    "%% 载入数据\n",
    "p=xlsread('粮食产量的8个影响因素.xls');\n",
    "t=xlsread('粮食总量和谷物、非谷物共3个量.xlsx');\n",
    "t1=1995:1:2009;\n",
    "t2=1995:1:2014;\n",
    "%训练数据\n",
    "P_train=p(1:15,:)';%训练输入数据\n",
    "T_train=t(1:15,:)';%训练输出数据\n",
    "\n",
    "%测试数据\n",
    "P_test=p(1:20,:)';%测试输入数据\n",
    "T_test=t(1:20,:)';%测试输出数据\n",
    "\n",
    "%数据归一化\n",
    "%输入数据归一化\n",
    "\n",
    "[Pn_train,inputps] = mapminmax(P_train,0,1);%原始训练数据归一化到0到1之间\n",
    "Pn_test = mapminmax('apply',P_test,inputps);%测试输入数据归一化到0到1之间\n",
    "% 输出数据归一化\n",
    "[Tn_train,outputps] = mapminmax(T_train,0,1);%原始输出数据归一化到0到1之间\n",
    "Tn_test = mapminmax('apply',T_test,outputps);%测试输出数据归一化到0到1之间\n",
    "\n",
    "%建立网络\n",
    "net=newff(minmax(Pn_train),[8,3],{'logsig','purelin'},'trainlm');\n",
    "\n",
    "net=train(net,Pn_train,Tn_train);\n",
    "%训练输入的预测结果\n",
    "Tn_output= sim(net,Pn_train);\n",
    "%仿真结果反归一化\n",
    "T_output = mapminmax('reverse',Tn_output,outputps);%反归一化的结果\n",
    "%\n",
    "figure\n",
    "plot(t1,T_train(1,:),'r-o','linewidth',2)\n",
    "hold on\n",
    "plot(t1,T_output(1,:),'b-diamond','linewidth',2)\n",
    "xlabel('训练样本')\n",
    "ylabel('粮食总量(万吨)')\n",
    "legend('实际值','BP预测值')\n",
    "figure\n",
    "plot(t1,T_train(2,:),'r-o','linewidth',2)\n",
    "hold on\n",
    "plot(t1,T_output(2,:),'b-diamond','linewidth',2)\n",
    "xlabel('训练样本')\n",
    "ylabel('谷物产量(万吨)')\n",
    "legend('实际值','BP预测值')\n",
    "figure\n",
    "plot(t1,T_train(3,:),'r-o','linewidth',2)\n",
    "hold on\n",
    "plot(t1,T_output(3,:),'b-diamond','linewidth',2)\n",
    "xlabel('训练样本')\n",
    "ylabel('非谷物产量(万吨)')\n",
    "legend('实际值','BP预测值')\n",
    "%%BP仿真\n",
    "Tn_sim = sim(net,Pn_test);%测试数据的仿真结果\n",
    "%仿真结果反归一化\n",
    "T_sim = mapminmax('reverse',Tn_sim,outputps);%仿真结果反归一化\n",
    "\n",
    "figure\n",
    "plot(t2,T_test(1,:),'r-o','linewidth',2)\n",
    "hold on\n",
    "plot(t2,T_sim(1,:),'b-diamond','linewidth',2)\n",
    "xlabel('测试样本')\n",
    "ylabel('粮食总量(万吨)')\n",
    "legend('实际值','BP预测值')\n",
    "box off\n",
    "figure\n",
    "plot(t2,T_test(2,:),'r-o','linewidth',2)\n",
    "hold on\n",
    "plot(t2,T_sim(2,:),'b-diamond','linewidth',2)\n",
    "xlabel('测试样本')\n",
    "ylabel('谷物产量(万吨)')\n",
    "legend('实际值','BP预测值')\n",
    "box off\n",
    "figure\n",
    "plot(t2,T_test(3,:),'r-o','linewidth',2)\n",
    "hold on\n",
    "plot(t2,T_sim(3,:),'b-diamond','linewidth',2)\n",
    "xlabel('测试样本')\n",
    "ylabel('非谷物产量(万吨)')\n",
    "legend('实际值','BP预测值')\n",
    "box off\n",
    "\n",
    "figure\n",
    "plot(t1,t_train(1,:),'r-o','linewidth',2)\n",
    "hold on\n",
    "plot(t1,train_out(1,:),'b-diamond','linewidth',2)\n",
    "xlabel('训练样本')\n",
    "ylabel('粮食总量(万吨)')\n",
    "legend('实际值','GRNN预测值')\n",
    "figure\n",
    "plot(t1,t_train(2,:),'r-o','linewidth',2)\n",
    "hold on\n",
    "plot(t1,train_out(2,:),'b-diamond','linewidth',2)\n",
    "xlabel('训练样本')\n",
    "ylabel('谷物产量(万吨)')\n",
    "legend('实际值','GRNN预测值')\n",
    "figure\n",
    "plot(t1,t_train(3,:),'r-o','linewidth',2)\n",
    "hold on\n",
    "plot(t1,train_out(3,:),'b-diamond','linewidth',2)\n",
    "xlabel('训练样本')\n",
    "ylabel('非谷物产量(万吨)')\n",
    "legend('实际值','GRNN预测值')\n",
    "%测试数据仿真\n",
    "test_out=sim(net,pn_test);\n",
    "test_out=postmnmx(test_out,mint,maxt);\n",
    "figure\n",
    "plot(t2,t_test(1,:),'r-o','linewidth',2)\n",
    "hold on\n",
    "plot(t2,test_out(1,:),'b-diamond','linewidth',2)\n",
    "xlabel('测试样本')\n",
    "ylabel('粮食总量(万吨)')\n",
    "legend('实际值','GRNN预测值')\n",
    "box off\n",
    "figure\n",
    "plot(t2,t_test(2,:),'r-o','linewidth',2)\n",
    "hold on\n",
    "plot(t2,test_out(2,:),'b-diamond','linewidth',2)\n",
    "xlabel('测试样本')\n",
    "ylabel('谷物产量(万吨)')\n",
    "legend('实际值','GRNN预测值')\n",
    "box off\n",
    "figure\n",
    "plot(t2,t_test(3,:),'r-o','linewidth',2)\n",
    "hold on\n",
    "plot(t2,test_out(3,:),'b-diamond','linewidth',2)\n",
    "xlabel('测试样本')\n",
    "ylabel('非谷物产量(万吨)')\n",
    "legend('实际值','GRNN预测值')\n",
    "box off\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "0bf4da5b",
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3 (ipykernel)",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.8.10"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 5
}
